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Abstract 


A  numerical,  finite  difference  computer  program  based  on  the  Navier- 
Stokes  equations  is  modified  to  give  avalanche  runout  distance,  velocity 
of  the  leading  edge  of  the  avalanche  and  depth  of  debris  in  the  runout 
zone.  The  program  requires  a  longitudinal  profile  of  the  avalanche  path, 
thickness  of  the  snow  in  the  starting  zone,  and  two  friction  coefficients. 
Kineomatic  viscosity  (or  coefficient  of  internal  friction)  is  expressed  as 
an  average  for  the  entire  avalanche  path.  Slope  inclinations  and 
coefficient  of  surface  friction  can  be  varied  for  every  10-  to  20-m 
increment  of  the  avalanche  path.  The  program  can  be  modified  to  allow 
for  variable  width  of  flow  and  for  snow  entrainment  during  flow.  In  test 
cases,  the  program  predicted  results  consistent  with  observations  for  a 
number  of  avalanche  events. 
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Numerical  Simulation  of 
Snow  Avalanche  Flow 

T.E.  Lang,  K.L.  Dawson,  and  M.  Martinelli,  Jr. 


Management  Implications 

Many  homes,  lodges,  and  other  structures  are  now 
being  built  in  steep,  snow- covered  areas  of  the 
western  U.  S.  Snow  avalanches  are  a  threat  at  many  of 
these  sites.  Land  managers,  county  commissioners, 
and  other  local  authorities  need  information  on  the 
location  of  avalanche  paths  and  the  frequency  and 
intensity  of  avalanches  to  safeguard  the  public  and 
still  permit  orderly  growth  of  mountain  communities. 

Potential  development  sites  should  be  examined 
and  avalanche  zoning  maps  prepared  before  land 
values  and  political  pressures  increase.  To  do  this,  a 
reliable  and  objective  way  is  needed  to  predict  the 
runout  distance,  impact  force,  and  return  interval  of 
avalanches  of  various  sizes  on  specific  avalanche 


paths.  From  these  data,  avalanche  areas  can  be 
mapped  on  the  basis  of  the  "50-year"  or  "100-year" 
avalanche  in  the  same  way  floodplains  are  mapped 
for  the  "100-year"  flood. 

This  paper  presents  a  new  technique  for  estimat- 
ing runout  distance  and  debris  depths  for  specific 
avalanche  paths  based  on  terrain  and  snow  condi- 
tions. Impact  forces  and  return  intervals  will  not  be 
discussed. 

Numerous  approaches  to  avalanche  dynamics 
have  been  undertaken  in  the  past,  and  new  or 
modified  approaches  wiU  undoubtedly  continue  to  be 
developed  to  give  better  solutions  to  the  complex 
problems  associated  with  avalanche  movement. 


Introduction 

The  prediction  of  runout  distance  is  an  important 
factor  in  establishing  avalanche  zoning  restrictions 
in  many  mountain  regions.  An  approximate, 
strongly  empirical,  hydrodynamic  method  for  pre- 
dicting runout  distance  developed  by  Voellmy  (1964) 
has  been  used  as  a  basis  for  avalanche  zoning  in 
Switzerland  for  several  years.  The  method  is  based 
on  uniform,  open-channel-flow  hydrodynamics, 
which  relates  flow  velocity  to  channel  slope  and 
dimensions  and  to  surface  roughness.  The  condition 
of  equilibrium  flow  restricts  the  application  of  this 
approach  to  slowly  varying  changes  in  the  govern- 
ing parameters.  The  two  basic  parameters  appearing 
in  the  velocity  equation  are  a  turbulence  coefficient 
and  a  surface  friction  coefficient.  These  parameters 
are  varied  to  take  into  account  surface  friction,  slope 
irregularity,  object  instrusions  (e.g.,  trees  and  rock 
outcrops),  flow  depth  changes,  and  other  phenomena 
of  the  flow  that  are  characteristic  of  avalanches.  In 
general,  avalanche  dynamists  agree  that  the 
successful  application  of  Voellmy's  equations  re- 
quires a  great  deal  of  professional  judgment  which 
can  only  be  gained  by  repeated  use  of  the  equations. 

Recognizing  the  transient  behavior  of  the  flow,  and 
the  need  to  separate  changes  in  path  geometry  from 
the  flow  parameters,  a  numerical  approach  appears 
to  be  a  better  way  to  represent  actual  avalanche 


motion  history  than  the  open  channel  flow  technique. 
The  basic  flow  parameters  can  then  be  varied  at  each 
increment  along  the  fluid  path  adding  flexibility  to 
the  model.  The  Los  Alamos  Scientific  Laboratory  has 
considerable  experience  in  transient  flow  using  finite 
differencing  of  mass  and  momentum  transfer  equa- 
tions (Amsden  and  Harlow  1970,  Welch  et  al  1966, 
Hirt,  Nichols  and  Romero  1975).  The  recently  re- 
ported simplified  fluid  dynamics  code,  SOLA,  (Hirt, 
Nichols  and  Romero  1975)  has  the  essential  features 
needed  to  model  avalanche  flow.  These  include 
capability  to  represent  a  free  surface  geometry 
(SOLASURF),  to  incorporate  nonzero  fluid  viscosi- 
ty, and  to  model  a  varying  lower  boundary  profile. 
Major  questions  that  are  not  answerable  from  cases 
considered  in  the  Los  Alamos  references  include 
whether  or  not  long- duration  flow,  such  as  in  an 
avalanche,  remains  stable,  and  whether  an  avalanche 
can  be  economically  modeled  and  analyzed.  This 
paper  examines  these  questions  and  modifies  the 
SOLA  code  to  permit  the  modeling  of  avalanche 
runout  distance. 

This  report  also  describes  in  detail  modifications 
to  the  SOLASURF  code  to  adapt  it  to  the  avalanche 
runout  problem.  The  modified  code  is  then  applied  to 
a  number  of  avalanche  cases  to  establish  bounds  on 
the  kinematic  (internal)  viscosity  and  the  surface 
friction  coefficient. 
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Description  of  the  Avalanche  Computer  Code 

Basic  Features  of  the  SOLASURF  Code 

The  governing  equations  of  motion  incorporated  in 
the  SOLASURF  code  in  finite  difference  form  are  the 
Navier- Stokes  equations: 
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In  equation  [1],  u  and  v  are  velocity  components;  g^ 
and  gy  are  gravity  components;  p,  is  the  density  of 
the  incompressible  fluid  and  is  combined  with  pres- 
sure, p,  to  become  the  actual  term  solved  for  in  the 
code  and  referred  to  as  "pressure;"  v  is  the  kine- 
matic viscosity  (or  internal  friction).  The  two  linear 
operators  of  equation  [1]  are 
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Additionally,  we  impose  the  condition  of  conserva- 
tion of  mass  for  an  incompressible  fluid  given  by 
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Consistent  with  the  usual  finite  difference  ap- 
proach, the  fluid  domain  is  divided  into  a  uniform 
grid  of  rectangular  cells,  surrounded  by  a  single-cell 
thickness  of  boundary  cells.  The  notations  associat- 
ed with  the  grid  layout  are  shown  in  figure  1 . 

In  the  code  IBAR^  and  JBAR  are  the  number  of 
flow  domain  cells  in  the  x  and  y  directions,  respec- 
tively, and  IMAX  and  JMAX  are  the  corresponding 
number  of  cells  including  the  boundary  cells.  Five  2- 
dimensional  arrays  with  the  dimensions  IMAX  by 
JMAX,  establish  the  major  memory  allocation.  The 
arrays  are  the  current  velocity  components  u  and  v, 
the  updated  velocity  components  u^^  and  v^,  and  the 
pressure,  p,  in  each  cell.  Limiting  the  size  of  these 
arrays  is  the  primary  method  of  controlling  running 
time  on  the  computer. 

'^Names  in  capital  letters  refer  to  parameters  used  in  the 
SOLASURF  code. 
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Figure  1  .—Grid  layout. 

The  algorithm  starts  with  an  initial  distribution  of 
fluid  in  the  domain  and  computes  a  pressure  and 
velocity  distribution  using  a  finite  difference  version 
of  equation  [1].  In  general,  the  computed  velocity 
components  do  not  satisfy  equation  [3],  so  that  small 
adjustments  of  cell  pressure  are  required  to  achieve  a 
zero-divergence  condition  to  a  prescribed  level  of 
accuracy  controlled  by  parameters  EPSI  and  DZRO 
(Hirt,  Nichols  and  Romero  1975).  The  initial  calcula- 
tion of  velocity  and  pressure  using  equation  [1]  is  the 
start  of  a  cycle  (CYCLE)  of  calculations,  and  each 
small  adjustment  of  parameters  in  one  sweep  of  the 
cells  is  an  iteration  (ITER). 

In  the  SOLASURF  code,  a  choice  between  four 
types  of  boundary  conditions  is  given  for  each 
boundary.  In  the  modified  version,  a  single  type  of 
boundary  condition  is  imposed  along  each  edge.  The 
left  wall,  considered  the  crown  face  of  the  avalanche, 
is  designated  a  rigid,  slip-free  boundary  in  the 
notation  of  Hirt,  Nichols  and  Romero  (1975).  The 
right  wall,  which  is  placed  beyond  the  runout  termi- 
nus, is  a  continuative  outflow  boundaiy,  as  is  the  top 
boundary.  The  bottom  boundary  is  modified  from 
the  rigid  no-slip  or  free-slip  options  given  in 
SOLASURF  to  model  different  friction  conditions 
that  might  be  encountered  in  an  avalanche  flow. 
Details  of  the  various  aspects  of  the  modified  pro- 
gram are  given  in  the  next  section. 

Output  from  SOLASURF  includes  two  types  of 
data,  depending  upon  user  preference.  The  more 
complete  output  at  the  end  of  an  integral  number  of 
cycles,  as  controlled  by  the  parameter  CWPRT,  is  a 
listing  of  the  velocity  and  pressure  in  each  cell.  At  the 
end  of  cycles  for  which  the  complete  printout  is  not 
requested,  a  single  line  is  printed  which  includes  the 
cycle  number  (CYCLE),  number  of  iterations  in  that 
cycle  (ITER),  time  increment  (DELT),  total  time  into 


2 


flow  (TIME),  fluid  volume  (FVOL),  maximum  slope- 
parallel  velocity  in  the  flow  (UMAX),  leading  edge 
slope-parallel  velocity  (UEDG),  and  cell  number  of 
the  leading  edge  of  the  avalanche  (LDEG). 

An  additional  input  parameter  is  TWFESf,  which  is 
the  time  estimate  for  the  avalanche  flow.  Normally, 
this  is  set  at  a  large  value  to  be  sure  the  complete 
avalanche  flow  is  represented  in  a  run.  However,  it 
may  also  be  used  to  terminate  computations,  in  which 
case  an  extended  printout  of  the  flow  distribution  is 
output. 

Remaining  input  parameters  needed  for  program 
execution  are  detailed  in  the  next  section,  since  they 
are  associated  directly  with  major  modifications 
made  to  the  code.  A  complete  description  and  an 
example  of  the  input  data  required  are  given  in 
Appendix  A. 

The  name  adapted  for  the  modified  version  of 
SOLASURF  for  the  simulation  of  avalanche  flow  and 
runout  distance  is  program  AVALNCH,  which  will  be 
used  in  subsequent  references  to  the  code. 

The  program  is  set  up  to  use  a  specific  set  of  units. 
All  linear  dimensions  are  expressed  in  meters,  and 
time  in  seconds.  Velocity  has  the  units  m/s;  the 
pressure/density  ratio  has  the  units  m"/s";  the  kine- 
matic viscosity  has  units  m^/s,  and  acceleration  due 
to  gravity  is  set  at  9.80  m/s^. 


one  cell.  An  improvement  on  the  grid  of  figure  2  is 
that  shown  in  figure  3.  The  primary  difficulty  is  the 
sensitivity  of  the  governing  equations  in  AVALNCH 
to  height  changes  between  successive  cells  along  the 
path.  This  is  not  remedied  by  the  grid  of  figure  3.  The 
final  grid  layout,  which  appears  most  adaptable  to 
the  avalanche  problem,  is  shown  in  figure  4,  and  is  the 
one  currently  programmed.  For  each  cell,  the  value  of 
AH  (fig.  3)  is  input,  from  which  gravity  components 
are  calculated  by  the  equations 
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Modifications  of  the  SOLASURF  Code 

There  are  several  limitations  to  developing  a 
numerical  model  for  avalanche  runout  prediction. 
First,  material  properties  must  be  expressed  in 
nominal  values,  since  point  values  of  the  properties 
as  a  function  of  slope,  position  and  time  are  unknown. 
Second,  the  calculations  should  be  low-cost,  because 
computer  and  manpower  resources  are  often  limited. 
Third,  representation  of  the  flow  should  be  simple 
and  not  require  complex  interpretation  of  parameter 
types  and  ranges.  These  three  constraints  have  been 
met  in  the  AVALNCH  code  with  the  imposition  of 
several  major  modifications  that  are  described 
below. 

If  a  grid  array,  based  upon  the  original  require- 
ments of  SOLASURF,  is  fit  to  an  avalanche  slope 
profile,  it  appears  as  shown  in  figure  2.  This  involves 
a  sizable  number  of  elements  in  the  y  direction,  since 
at  no  location  can  the  slope  exceed  the  ratio  S  J8^. 
With  this  restriction,  8^,  is  of  the  same  order  as  5^  (10 
to  30  m)  for  slopes  of  from  30  to  50  degrees,  so  that  an 
avalanche  2  m  high  occupies  only  a  small  fraction  of 


Figure  2.— Avalanche  profile  and  grid. 


Figure  3.— Improved  grid  layout 
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Figure  4.— Equivalent  grid  layout  for  avalanche  flow. 
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Angle  0  is  the  average  angle  of  the  slope  at  the  cell 
being  specified.  With  this  arrangement,  Sy  can  be 
selected  independent  from  Sx,  where  Sx,  is  the  slope- 
parallel  increment  along  the  avalanche  path.  For 
local  instrusions,  such  as  a  level  road  bed,  we  set 
AH=0  for  the  cells  in  the  road  cut.  Li  some  cases  a 
road  cut  may  also  cause  an  increase  in  the  airborne 
component  of  a  flowing  avalanche  by  deflecting  flow 
upward  and  entraining  air.  This  can  be  represented 
in  program  AVALNCH  by  reducing  the  values  of  v 
and  f  (the  coefficient  of  friction  at  the  lower  boun- 
dary) on  the  basis  of  available  physical  data. 

The  complex  flow  associated  with  a  vertical  drop  is 
not  representable  in  the  present  formulation.  The  g^ 
and  g  values  enter  the  governing  equations  in 
AVAL^vJCH  at  the  stage  of  preliminary  calculation  of 
the  velocity  and  pressure  distribution  at  the  start  of 
each  cycle.  The  ratio  Sy/^^  enters  later  when  bounda- 
ry conditions  are  imposed,  and  when  the  updated 
flow  height  is  computed.  By  numerical  experimenta- 
tion, it  is  determined  that,  whereas  the  limit  of  8 JS^ 
<  1.0  is  a  sensitive  parameter  to  numerical  stability, 
values  of  gy/g^  1-0  are  readily  handled.  The 
difference  is  that  rapid  variation  in  velocity  between 
cells  is  permissible,  but  rapid  variation  in  boundary 
height  is  not.  The  values  of  AH(HN(I))  for  each  cell 
are  input  as  a  one- dimensional  array  following  spe- 
cification of  the  snow  height,  H(I),  also  a  one-dimen- 
sional array. 

Velocity  at  the  lower  boundary,  u„,  is  determined 
by  linear  interpolation  between  the  corresponding 
velocities  in  the  boundary  cell,  u-j^,  and  the  lowest 
flow-domain  cell  U2,  (fig.  5)  in  which  velocity  u^  is 
expressed  in  terms  of  velocity  U2.  We  generalize  this 
relationship  from  that  of  the  SOLASURF  code  by 
specifying  tnat 


Up(1  -f) 
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Ui  =U2(1  -2f) 


[5] 


which  is  equivalent  to  a  surface  velocity  given  by 


Forf  =  0,  u-|^=  U2,  and  the  surface  is  slip-free,  and  for 
f  =  1,  u-j^  =  — U2,  so  that  u  g  =  0  a  no-slip  condition.  For 
intermediate  values  of  f,  the  boundary  is  a  partial  slip 
wall.  Since  f  reflects  friction  at  the  boundary,  we  refer 
to  this  as  the  "friction  coefficient,"  which  is  one  of 
two  basic  parameters  that  must  be  sized  for  the 
avalanche  runout  problem. 

Two  options  for  input  of  f  are  given  in  AVALNCH. 
If  f  is  constant  over  the  entire  path,  then  this  value  is 
input  by  specifying  parameter,  FRK.  If  f  is  a  variable 
over  the  path,  then  a  one-dimensional  array,  FRC(I), 
is  specified,  and  FRK=  0.0.  If  the  friction  is  actually 
to  equal  zero,  it  must  be  input  by  the  array.  If ,  in  a 
particular  application,  values  of  f  >  1.0  are  used,  this 
can  be  interpreted  as  a  penetration  of  the  surface 
roughness  into  the  flow.  The  surface  roughness 
where  u  =  0.0  is  defined  by 


r  = 


f-1 
2f 
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as  denoted  in  figure  5.  There  is  an  artificial  aspect  to 
this  interpretation  in  that  the  quantity  of  flowing  fluid 
is  not  reduced  correspondingly  due  to  r  >  0.  Thus, 
the  analogy  that  perhaps  most  closely  represents  the 
case  f  >  1.0  is  sparsely  spaced  rock  protrusions  or 
tree  trunks  that  slow  the  flow  significantly,  but  do  not 
trap  large  fractions  of  it.  The  variation  in  r  for 
different  values  of  f,  as  predicted  by  equation  [7],  is 
shown  in  figure  6. 

The  other  parameter  in  AVALNCH  requiring 
specification  is  the  kinematic  viscosity,  v  (NU).  It  is 
input  as  a  single  constant  over  the  entire  flow  path, 
since  information  on  possible  variation  of  this  para- 
meter is  generally  lacking.  When  we  learn  enough  to 
be  able  to  specify  changes  in  this  variable  down  the 
path,  program  AVALNCH,  as  written,  can  easily 
accommodate  this. 
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Figure  5.— Effect  of  friction  on  boundary  flow. 
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Figure  6.— Obstacle  protrusion  factor,  r,  versus  friction 
coefficient  f>1.0. 

From  case  studies  it  has  been  determined  that  the 
internal  circulation  in  a  flow  is  a  function  of  the 
number  of  cells  used  in  the  y  direction,  for  a  fixed 
number  of  cells  in  the  x  direction.  Since  values  for  f 
and  V  must  be  established  from  case  studies  rather 
than  from  controlled  experiments,  their  absolute 
values  are  not  important.  Thus,  to  standardize  the 
cell  geometry  for  all  avalanche  runs,  a  single  cell 
height  of  the  same  dimension  as  the  maximum 
thickness^  of  the  fracture  face  is  used.  Flow  height^ 
as  given  by  program  AVALNCH,  is  the  height  of  the 
"core  material"  and  not  of  the  snowdust  cloud  that 
accompanies  mixed  flow  and  powder  avalanches. 
Actual  field  observations  of  "core  height'"  are  very 
difficult  because,  in  most  cases,  only  the  dust  cloud 
can  be  seen  or  photographed.  A  second  flow  domain 
cell  is  assumed  above  the  avalanche  cell,  so  that 
JBAR  =  2  for  all  cases.  Based  upon  this  selection  of  a 
flow  domain  geometry,  f  and  u  are  sized  from  a 
number  of  modeled  avalanche  runs. 

Besides  numerical  stability  problems  associated 
with  the  variable  boundary  option  in  SOLASURF, 
another  difficulty  is  flow  of  mass  through  the  bounda- 
ry. Initial  avalanche  runs  showed  gradual  disappear- 
ance of  mass  until,  at  some  point  downstream,  the 
avalanche  left  the  flow  domain  completely.  This 
problem  is  reduced  by  using  the  standardized  grid  of 
figure  4,  and  is  eliminated  by  inserting  a  multiplica- 


^7^76  term  "thickness"  is  used  to  designate  the  dimension  of  a 
snow  layer  perpendicular  to  the  sliding  surface.  When  speak- 
ing of  the  flowing  snow,  however,  the  term  "flow height"  Is  used 
to  describe  the  dimension  perpendicular  to  the  slope  to  avoid 
using  the  term  "flow  thickness,"  which  could  be  misunder- 
stood as  relating  to  the  viscosity  of  the  flowing  snow. 


five  factor  equal  to  the  ratio  of  original  volume  to 
current  volume  in  the  height  update  section  (4000 
section)  of  the  code  (Appendix  B).  This  ratio  is 
close  to  unity,  so  that  the  correction  on  each  cycle  is 
small,  but  sufficient  to  eliminate  the  monotonic 
divergence  obsei^^ed  without  the  correction. 

One  difficulty  of  the  circulating  flow  condition  is  in 
obtaining  an  accurate  value  for  the  average  advance 
or"group"  velocity  of  the  fluid.  (A"group"  velocity  is 
needed  rather  than  a  particle  velocity.)  The  particle 
velocity,  which  reflects  circulation  as  well  as  transla- 
tion, listed  in  the  expanded  output  of  the  flow  profile 
is  the  maximum  value  of  pai'ticle  velocity  (UAIAX) 
and  is  used  to  predict  the  next  time  increment 
(BELT)  in  AVALNCH.  To  obtain  the  "group"'  veloc- 
ity in  a  mathematically  rigorous  manner  requires  a 
lengthy  calculation.  Instead,  an  approximate  calcula- 
tion utilizing  the  ceU  number  of  the  leading  (LDEG) 
and  trailing  (KTEG)  edges  of  the  avalanche  is  used. 
These  cells  are  designated  as  the  first  occurring  from 
opposite  directions  in  which  the  flow  height  is  iSc  (or 
less)  of  the  original  avalanche  height.  Once  the 
leading  edge  is  located,  the  time  change  for  advance 
of  this  edge  to  the  next  ceU  divided  into  the  cell  length 
is  one  estimate  of  the  avalanche  slope-pai'allel  veloci- 
ty (UEDG).  Since  fraction  distance  into  the  cells  is 
not  known,  a  nominal  velocity  estimate  must  be 
based  upon  averaging  the  values  of  UEDG  over 
several  (three  to  five)  cells.  Furthermore,  UEDG  can 
not  be  computed  until  the  leading  edge  has  left  the 
cell,  so  that  the  listed  values  of  UEDG  must  actually 
be  associated  with  the  cell  previous  to  the  current 
leading  edge  cell  (LDEG). 

A  condition  observed  in  early  runs  on  avalanches 
was  a  tendency  of  the  flow  to  reach  a  near  steady- 
state  shallow  flow,  not  unlike  equilibrium  open  chan- 
nel flow,  and  to  run  out  for  considerable  distances. 
Since  this  is  not  observed  in  most  avalanches,  a 
velocity  dependent  surface  friction  law  is  incorporat- 
ed in  the  code  (6000  section)  in  which  the  friction 
coefficient  is  increased  exponentially  with  decrease 
in  speed.  The  relationship  established  between  f  and 
u  is: 


f  =fo  (1  +20e 


-1.25  u\ 


[8] 


where  f^  is  the  initial  value  selected  for  f.  This  surface 
friction- velocity  relationship  is  wTitten  into  the  pro- 
gram to  be  operable  anytime  after  the  primary 
acceleration  phase  of  the  avalanche  flow.  A  plot  of 
equation  [8]  showing  the  variation  in  f/f^  versus  u  is 
given  in  figure  7.  Although  an  actual  quantitative 
physical  basis  for  this  phenomenon  is  not  known,  this 
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Figure  7.— Surface  friction  coefficient  ratio  as  a  function  of 
flow  velocity  for  "Fast-Stop"  option  in  AVALNCH. 

"fast-stop"  option  is  a  better  representation  for  most 
avalanche  flows  than  what  is  observed  without  it. 

Apart  from  user  program  limits,  there  are  five 
default  conditions  within  AVALNCH  which  can 
terminate  calculations.  They  are: 

1 .  Number  of  iterations  (ITER)  within  a  cycle 
reaches  500. 

2.  Duration  of  avalanche  flow  exceeds  time 
specified  (TWFIN). 

3.  Leading  edge  flow  velocity  (UEDG)  drops 
to  less  than  5%  of  maximum  flow  velocity 
achieved  during  the  flow. 

4.  Leading  edge  of  flow  reaches  lower  limit  of 
grid. 

5 .  Leading  edge  of  flow  does  not  advance  to  a 
new  cell  in  50  cycles. 

Parameter  overflow  occurs  for  condition  1.  For  con- 
ditions 2  through  5,  the  flow  distribution  is  printed 
and  an  end  statement  is  given  prior  to  termination. 

These  are  the  major  modifications  incorporated  in 
program  AVALNCH  Other  minor  modifications  have 
been  made  that  can  be  identified  by  comparing  the 
listing  of  AVALNCH  (Appendix  B)  to  that  of 
SOLASURF  (Hirt,  Nichols  and  Romero  1975). 


Characteristics  of  Fluid  Modeling  Equations 

Having  reviewed  the  various  features  that  have 
been  incorporated  in  program  AVALNCH,  we  con- 
sider the  physical  conditions  that  are  modeled  in  the 


computer  formulation.  The  Navier- Stokes  equations 
are  equations  of  motion  for  laminar  fluid  flow.  If 
internal  viscosity  is  greater  than  zero,  then,  the  flow  is 
also  considered  to  be  rotational.  In  laminar,  rotation- 
al flow  fluid  particles  tend  to  follow  curved  paths 
along  streamlines,  and  circulation  patterns  are 
evident.  In  contrast,  particle  motion  in  turbulent  flow 
is  completely  random,  and  streamlines  or  flow  path 
lines  are  not  evident.  From  the  study  of  numerous 
avalanche  films,  our  conclusion  is  that  for  flowing  and 
mixed  motion  avalanches  the  flow  is  often  rotational 
and  laminar.  See  Perla  and  MartinelH  (1976)  for 
description  of  avalanche  motion. 

At  the  extremes  of  avalanche  motion  are  the  cases 
of  sliding  rigid-blocks  (sliding  motion)  for  which 
internal  viscosity  is  large,  and  turbulent  flow  (powder 
avalanches)  for  which  viscosity  is  low.  The  low 
viscosity  of  turbulent  flow  is  attributable  to  a  large 
air-fraction  in  the  flow  and  small  snow  particles, 
which  produces  the  snow  dust  clouds  surrounding 
the  larger  mixed  motion  avalanches.  However,  even 
in  this  case,  there  is  a  core  of  denser  snow  within  the 
snow  dust  cloud  that  acts  as  a  piston  to  drive  or 
propel  the  airborne  component.  For  this  type  of 
avalanche,  program  AVALNCH  can  be  used,  but  it 
models  only  the  flow  of  the  core  material,  not  the 
snow  dust  cloud.  Once  the  snow  is  in  motion, 
subsequent  flow  height  predicted  by  the  program  is 
the  height  of  the  core  material,  not  that  of  the  snow 
dust  cloud. 

Most  field  observations  of  the  flow  height  of  this 
type  of  avalanche  are  to  the  top  of  the  snow  dust 
cloud,  which  gives  no  information  on  the  type  of  flow 
that  the  core  exhibits.  However,  based  upon  evi- 
dence from  debris,  the  assumption  of  laminar  rota- 
tional flow  of  the  core  material  appears  reasonable  in 
many  cases. 

Program  AVALNCH  traces  the  snow  flow  starting 
from  the  initial,  at-rest  conditions  of  the  release  slab 
in  the  starting  zone.  Specification  of  the  distribution 
of  snow  slab  in  the  starting  zone  is,  thus,  a  necessary 
part  of  the  input.  Any  slab  distibution  may  be  input, 
the  only  restrictions  are  the  following: 

1.  Only  one  height  per  slope-parallel  cell  may 
be  specified;  this  height  is  assumed  to  be 
constant  over  the  entire  cell  and  is  mea- 
sured perpendicular  to  the  slope. 

2.  Cell  height  is  then  set  equal  to  the  maxi- 
mum slab  thickness  specified  in  the  start- 
ing zone  cells. 


6 


I 


Snow  entrainment  as  flow  propagates  down  the 
slope  is  common  in  avalanche  flow.  To  model  this 
phenomenon,  it  is  not  accurate  to  specify  snow 
thickness  in  downslope  cells  in  the  same  way  that 
starting  zone  material  is  specified.  However,  the 
modeling  of  snow  entrainment  is  possible  with  pro- 
gram AVALNCH  and  is  considered  in  detail  in  the 
section  headed  "Leading  Edge  Snow  Entrainment 
Modeling." 


Avalanche  Models,  and  Evaluations 
Using  Program  AVALNCH 

Intent  of  Case  Studies 

This  part  of  the  paper  outlines  a  number  of 
specific,  recorded  avalanche  incidents  to  establish 
ranges  for  kinematic  viscosity,  v,  and  surface  fric- 
tion, f.  Ideally,  these  ranges  would  be  established 
using  detailed  information  from  many  avalanches 
based  on  the  following  factors:  fracture  thickness, 
slope-parallel  cross  section  of  the  initially  released 
slab,  type  of  surface  layer  of  snow  on  slope,  velocity 
of  flow,  actual  runout  distance,  and  distribution  of 
terminal  debris.  Since  such  detailed  data  are  seldom 
available,  we  used  case  studies  to  size  f  and  v  to 
model  overall  runout  response.  Such  a  general  ap- 
proach is  thought  to  have  greater  utility  than  one 
requiring  highly  accurate  point  variations  in  para- 
meters because  of  the  difficulty  in  obtaining  such 
data  for  a  wide  range  of  snow  and  slope  conditions. 

Any,  and  ideally  all  of  the  above  parameters,  if 
measured,  could  be  used  to  better  quantify  the 
values  of  v  and  f.  However,  the  tremendous  varia- 
tion in  snow  and  slope  conditions  that  occur  through- 
out a  winter  season  preclude  the  utility  of  a  specific 
slope  analysis  taking  into  account  point  variations  in 
parameters.  In  fact  a  program  requiring  this  accuracy 
would  not  be  of  general  applicability  for  runout 
prediction,  since  specification  of  parameters  is  sel- 
dom accurate  except  under  very  specialized  con- 
ditions. 

Ironton  Park  Avalanche  Path 

Location:  San  Juan  Mountains,  Red  Moun- 

tain Pass  North;  southeast  slope  of 
Hayden  Mountain,  Colorado. 

^The  term  "runout  distance"  is  used  to  indicate  ttie  horizontal 
distance  from  the  start  of  the  runout  zone  to  the  most  distant 
end  of  the  avalanche  debris.  In  this  paper,  the  term  "travel 
distance"  is  used  to  mean  the  horizontal  distance  from  the 
fracture  face  in  the  starting  zone  to  the  most  distant  end  of  the 
avalanche  debris. 


Starting  Zone: 


Track: 


Runout  Zone: 


Avalanches: 


Profile: 

Length  (m) 

0 
70 

174 
265 
326 
418 

509 
814 

997 

1,073 


Wide,  uniform  slope  below  Half 
Moon  Basin,  below  timberline;  ele- 
vation 3,230  to  3,450  m. 

The  entire  slope;  vertical  drop  350 
m;  length  1,000  m. 

Level  bottom  of  Ironton  Park  and 
frozen  lake. 

Large  dry  snow  avalanche  in  1958 
crossed  the  road  and  left  a  debris 
pile  approximately  1.5  m  high  on 
the  road.  No  information  is  known 
on  the  starting  zone  fracture  thick- 
ness or  extent,  nor  on  the  type  of 
snow  constituting  the  flow. 

Ironton  Quadrangle,  Colorado  (7.5'). 
Elevation  (m) 


354 
293 
232 
171 
110 
49 
0 
0 


0 


49 


Location  of 
road  (Cell  93) 
Opposite  side 
of  flat 


The  profile  constructed  from  these  data  points  is 
shown  in  figure  8.  The  data  needed  from  this  profile  is 
obtained  by  first  marking  off  the  increment  distance 
8^  along  the  profile.  This  is  selected  to  be  10  m  slope 
distance  in  this  case  in  order  to  describe  the  profile 
accurately  and  still  keep  the  number  of  cells  to  less 
than  the  limit  of  200  for  AVALNCH.  The  change  in 
elevation  in  each  cell  is  taken  directly  from  figure  8. 
For  the  present  example  only,  the  values  of  the  eleva- 
tion and  change  in  elevation  for  each  cell  are  listed  in 
table  1.  It  should  be  noted  that  positive  values  of  AH 
correspond  to  decreases  in  elevation  across  cells,  and 
negative  values  of  AH  correspond  to  increases  in 
elevation  across  cells  (adverse  slope).  The  values 
shown  in  table  1  for  AH,  carried  out  to  tenths  of  a 
meter,  were  obtained  by  interpolation,  and  reflect  an 
accuracy  not  necessary  for  meaningful  avalanche 
prediction.  In  subsequent  runs,  whole  number  values 
for  AH  were  used  successfully,  with  much  less  time 
spent  in  data  preparation.  However,  gross  round-off 
of  AH  to  the  point  of  masking  5-  to  20-m  length 


perturbations  in  slope  geometry  is  not  recommended, 
as  this  order-of-magnitude  perturbation  can  have  a 
significant  effect  on  the  flow.  Some  of  these  condi- 
tions will  be  pointed  out  in  the  following  avalanche 
cases. 

In  the  case  of  the  Ironton  Park  slope,  little  infor- 
mation is  known  about  specific  avalanches.  Since  the 
track  is  an  open  uniform  slope,  this  example  was  used 
to  examine  the  effect  of  varying  the  values  of  the 
kinematic  viscosity,  v,  and  surface  friction,  f.  A 
summary  of  the  values  of  v  and  f  used  for  different 
travel  distances  is  given  in  table  2.  For  most  of  the 
runs  the  initial  slab  thickness  is  set  at  2  m  for  a 
downslope  distance  of  50  m,  then  a  linear  decrease  in 
thickness  out  to  90  m,  which  is  considered  a  large 
initial  release.  Other  runs  reported  in  table  2  for 
initial  slab  thickness  other  than  the  2  m  are  vertically 
scaled  modifications  of  the  2-m  slab  distribution. 

The  avalanche,  which  occurred  in  late  December 
1958  and  deposited  1.5  m  of  snow  on  the  road,  is 
considered  rare.  Unfortunately,  fracture  line  infor- 
mation is  lacking.  Occurrence  of  the  avalanche  in  late 
December  suggests  moderate  midwinter  snow.  On 
the  basis  of  an  initial  slab  thickness  of  2  m,  computer 
runs  that  characterize  avalanche  termination  in  the 
vicinity  of  the  road  are  those  for  which  values  of  v 
and  f  sum  to  1.0  to  1.1.  For  example,  for  v  =  f  =  0.5, 
the  avalanche  leading  edge  is  60  m  beyond  the  road, 
and  the  debris  on  the  road  is  nominally  1.9  m  deep 
(Appendix  C).  For  v  =  0.4,  f  =  0.7,  the  leading  edge 
of  the  avalanche  reaches  the  road,  and  for  v  =  0.5,  f 


=  0.6,  the  leading  edge  is  20  m  short  of  the  road.  It  is 
apparent  that  by  modifying  v  and  f  by  relatively 
small  amounts  a  wide  variety  of  snow  conditions  can 
be  modeled.  We  further  note  from  the  results  in  table 
2  the  near  correspondence  of  results  when  the  values 
of  V  and  f  are  exchanged.  That  is,  for  v  =  0.4andf= 
0.7,  travel  distance  is  approximately  the  same  as  for 
v=0.7  andf=0.4.  Similar  results  are  obtained  for  all 
other  cases  where  v  is  not  equal  to  f  reported  in  table 
2.  This  indicates  the  nearly  equal  influence  of  these 
two  parameters  in  determining  travel  distance  when 
using  the  units  of  measure  selected  for  these  studies. 
Thus,  in  the  present  case,  we  consider  u  =  f  -  0.5  as 
representative  of  the  flow,  and  adapt  the  practice  of 
selecting   f  =  f ,  except  when  physical  conditions 
clearly  indicate  different  values  for  v  and  f  are 
warranted. 

The  computer  listing  of  the  input  and  output  data 
for  the  case  v  =  f  =  0.5  and  a  slab  thickness  of  2  m  is 
shown  in  Appendix  C.  The  listing  includes  the  eight 
discrete  input  parameters,  followed  by  the  flow 
height,  elevation  change  per  cell,  and  the  boundary 
friction  coefficients,  for  the  purpose  of  data  verifica- 
tion. Finally,  the  gravity  components,  g^  and  g^,  are 
printed  for  reference  purposes.  Any  adverse  grade  in 
the  slope  profile  shows  up  with  negative  values  in  the 
array  for  g  .  The  values  of  g  should  always  be 
negative.  tIiIs  completes  printout  of  all  input  data 
needed  for  AVALNCH. 

Outputs  are  the  single  line  printouts  of  the  ava- 
lanche flow  status  given  after  each  cycle  of  calcula- 
tion. By  scanning  the  column  of  numbers  labeled 
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Figure  8.— Longitudinal  profile  of  Ironton  Park  avalanche  path. 
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Table  1.  —  Ironton  Park  profile  elevation  data 


Elevation  AH 

Elevation  AH 

Elevation  AH 

Elevation  AH 

Elevation  AH 

Station 

(m\ 
\i"  1 

Station 

(m) 

Station 

Station 

station 

0 

357.0 

23 

226.2 

4.6 

46 

81 .4 

6.0 

69 

2.0 

0.9 

92 

0.0 

0.0 

1 

350.5 

6.5 

24 

221.4 

4.8 

47 

75.4 

6.0 

70 

1 .3 

0.7 

93 

0.0 

0.0 

2 

344.3 

6.2 

25 

216.4 

5.0 

48 

69.8 

5.6 

71 

0.8 

0.5 

94 

0.0 

0.0 

3 

338.1 

6.2 

26 

211.2 

5.2 

49 

64.6 

5.2 

72 

0.4 

0.4 

95 

0.0 

0.0 

4 

331.9 

6.2 

27 

205.8 

5.4 

50 

59.4 

5.2 

73 

0.1 

0.3 

96 

0.0 

0.0 

5 

325.4 

6.5 

28 

200.2 

5.6 

51 

54.4 

5.0 

74 

0.0 

0.1 

97 

0.0 

0.0 

6 

319.2 

6.2 

29 

194.4 

5.8 

52 

49.8 

4.6 

75 

0.0 

0.0 

98 

0.0 

0.0 

7 

313.0 

6.2 

30 

188.2 

6.2 

53 

45.6 

4.2 

76 

0.0 

0.0 

99 

0.0 

0.0 

8 

306.8 

6.2 

31 

182.0 

6.2 

54 

41 .4 

4.2 

77 

0.0 

0.0 

100 

0.0 

0.0 

9 

300.3 

6.5 

32 

175.8 

6.2 

55 

37.4 

4.0 

78 

0.0 

0.0 

101 

0.0 

0.0 

10 

293.8 

6.5 

33 

169.6 

6.2 

56 

33.4 

4.0 

79 

0.0 

0.0 

102 

0.0 

0.0 

11 

287.6 

6.2 

34 

163.1 

6.5 

57 

29.8 

3.6 

80 

0.0 

0.0 

103 

0.0 

0.0 

12 

281.8 

5.8 

35 

156.3 

6.8 

58 

26.2 

3.6 

81 

0.0 

0.0 

104 

0.0 

0.0 

13 

275.8 

6.0 

36 

148.9 

7.4 

59 

22.8 

3.4 

82 

0.0 

0.0 

105 

0.0 

0.0 

14 

270.0 

5.8 

37 

141 .7 

7.2 

60 

19.4 

3.4 

83 

0.0 

0.0 

106 

0.0 

0.0 

15 

264.4 

5.6 

38 

134.5 

7.2 

61 

16.4 

3.0 

84 

0.0 

0.0 

107 

1 .0 

-1 .0 

16 

5  fi 

39 

1?7  4 

7  1 

62 

13  fi 

85 

0  0 

0  0 

108 

?  0 

-1  0 

17 

253.6 

5.2 

40 

120.4 

7.0 

63 

11.1 

2.5 

86 

0.0 

0.0 

109 

4.0 

-2.0 

18 

249.0 

4.6 

41 

113.4 

7.0 

64 

8.9 

2.2 

87 

0.0 

0.0 

110 

7.0 

-3.0 

19 

244.8 

4.2 

42 

106.4 

7.0 

65 

7.1 

1.8 

88 

0.0 

0.0 

20 

240.6 

4.2 

43 

99.8 

6.6 

66 

5.5 

1.6 

89 

0.0 

0.0 

21 

236.0 

4.6 

44 

93.6 

6.2 

67 

4.1 

1.4 

90 

0.0 

0.0 

22 

230.8 

5.2 

45 

87.4 

6.2 

68 

2.9 

1.2 

91 

0.0 

0.0 

Table  2.  —  Ironton  Park  avalanche  travel  distance  for  different  values  of  viscosity  and  surface  friction  coefficients 


Run 
No. 

Viscosity 
(mVs) 

Friction 
coeff. 

Initial 
slab 
thickness 
(m) 

Travel 
Distance^ 
(m) 

U  max  Cell 
(m/s)  No. 

Comments 

1 

0.3 

0.3 

2.0 

63/71 

Avalanche  flows  out  of  grid  at  25  m/s. 

2 

0.3 

0.5 

2.0 

54/65 

Avalanche  flows  out  of  grid  at  16  m/s. 

3 

0.4 

0.4 

2.0 

53/64 

Avalanche  flows  out  of  grid  at  1 5  m/s. 

4 

0.4 

0.5 

2.0 

47/60 

Avalanche  flows  out  of  grid  at  8  m/s. 

5 

0.4 

0.7 

2.0 

820 

39/58 

Avalanche  stops  on  road.  Max.  debris  depth:  2.4  m. 

6 

0.45 

0.45 

2.0 

47/62 

Avalanche  flows  out  of  grid  at  8  m/s. 

7 

0.45 

0.5 

2.0 

44/59 

Avalanche  flows  out  of  grid  at  2.5  m/s. 

8 

0.5 

0.5 

2.0 

880 

42/58 

Avalanche  stops  60  m  beyond  road.  Max.  debris 

depth  2.2  m. 

9 

0.5 

0.6 

2.0 

790 

Avalanche  stops  20  m  in  front  of  road. 

10 

0.5 

0.7 

2.0 

740 

33/55 

Avalanche  stops  70  m  in  front  of  road.  Max.  debris 

depth  2.1  m. 

11 

0.6 

0.5 

2.0 

790 

37/56 

Avalanche  stops  20  m  in  front  of  road.  Max.  debris 

depth  2.2  m. 

12 

0.6 

0.6 

2.0 

730 

33/55 

Avalanche  stops  80  m  in  front  of  road.  Max.  debris 

depth  2.2  m. 

13 

0.6 

0.7 

2.0 

700 

33/55 

Avalanche  stops  110  m  in  front  of  road.  Max.  debris 

depth  2.0  m. 

14 

0.7 

0.4 

2.0 

820 

39/63 

Avalanche  stepson  road.  Max.  debris  depth  2.0  m. 

15 

0.7 

0.5 

2.0 

740 

39/55 

Avalanche  stops  70  m  in  front  of  road.  Max.  debris 

depth  2.7  m. 

16 

0.7 

0.6 

2.0 

710 

29/54 

Avalanche  stops  100  m  in  front  of  road.  Max.  debris 

depth  2.0  m. 

17 

0.4 

0.4 

1.5 

820 

39/63 

Avalanche  stops  on  road.  Max.  debris  depth  1 .5  m. 

18 

0.45 

0.45 

0.5 

720 

33/56 

Avalanche  stops  90  m  in  front  of  road.  Max.  debris 

depth  2.4  m. 

19 

0.5 

0.5 

1.5 

690 

20/51 

Avalanche  stops  120  m  in  front  of  road. 

20 

0.5 

0.5 

1.0 

640 

14/50 

Avalanche  stops  50  m  onto  flat.  Max.  debris  depth 
1.4m. 

21 

0.5 

0.5 

0.75 

160 

After  initial  transient,  avalanche  decelerates 

continuously  and  stops  on  a  bench  of  the  slope. 

1  Travel  distance  is  the  horizontal  distance  from  the  crown  of  the  avalanche  to  the  leading  edge  of  the  terminal  debris 
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UEDG  and  averaging  over  several  cells,  we  see  that 
the  maximum  velocity  achieved  by  the  avalanche  is 
about  42  m/s  in  the  vicinity  of  cell  64.  This 
corresponds  to  a  distance  approximately  530  m  (fig. 
8),  which  is  the  start  of  the  runout  zone.  The 
avalanche  stops  in  cell  99,  one  cell  short  of  the  grid 
limit.  The  final  listing  of  data  is  the  extended  output. 
By  scanning  the  column  labeled  H  (Appendix  C),  we 
find  the  final  position  of  debris  to  extend  from  cell  91 
(about  20  m  in  front  of  road),  where  the  depth  of 
debris  is  0.8  m,  to  cell  99  where  the  debris  depth  is  0.02 
m  (Appendix  C).  Although  velocity  components 
persist  in  columns  for  u  and  v,  these  must  be 
considered  components  of  the  circulation  and  are 
assumed  to  be  negligible.  Values  for  velocities  and 
pressure  outside  the  range  of  the  flow  domain  can  be 
ignored. 

Holding  V  and  f  constant  at  0.5,  the  thickness  of 
the  initial  slab  is  varied  from  0.5  to  2.0  m,  and  travel 
distance  is  noted  as  shown  in  figure  9.  These  results 
show  a  nonlinear  dependence  of  runout  distance  with 
slab  height  with  the  "fast- stop"  option  having  a 
significant  effect  on  flow  of  the  shallower  avalanches. 

An  avalanche  with  a  2-m  average  slab  thickness  is 
considered  unusual.  What  is  perhaps  more  realistic  is 
a  thickness  of  1.5  m,  for  which  u  =  f  =  0.4  in  order 
that  the  debris  reach  the  road  (table  2,  run  17).  Thus, 
we  conclude  that  sensitive  parameters  that  influence 
avalanche  travel  distance  are  the  viscosity  and  fric- 
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tion  coefficients  and  initial  slab  thickness.  Correla- 
tion between  these  parameters  is  expected  as  addi- 
tional avalanche  occurrences  are  evaluated. 

The  principal  conclusions  we  obtain  from  this 
study  of  the  Ironton  Park  profile  are: 

1.  Values  of  v  and  f  in  the  range  0.4  to  0.5 
represent  nominal  values  for  early  winter,  dry 
hardpack  snow  in  thick- slab  avalanches. 

2.  Cell  lengths  of  the  order  of  tens-of-meters  in 
the  slope-pai'allel  direction,  coupled  with  cell 
heights  of  1  to  2  m  in  the  slope-normal 
direction,  can  be  selected  to  reasonably 
model  the  avalanche  flow. 

3.  Parameters  that  are  identified  to  strongly 
influence  avalanche  travel  distance  are  the 
kinematic  viscosity,  friction  coefficient,  and 
initial  slab  thickness.  Ranges  of  values  for 
these  parameters  should  be  established  as 
additional  case  studies  of  specific  avalanche 
paths  are  modeled. 


Pallavicini  Avalanche  Path 


Location: 


Starting  Zone: 


Track: 


Runout  Zone: 


Avalanches: 


Profile: 


Figure  9.— Avalanche  travel  distances  versus  initial  slab 
height  for  the  Ironton  Park  path,  holding  viscosity  and 
friction  coefficients  constant. 


Front  range  of  Colorado  Rockies; 
north  slope  of  summit  point  3,700  m 
southwest  of  Arapaho  Basin  ski 
area. 

Bowl-shaped  depression  in  the 
north  facing  slope;  above  and  below 
timberline;  elevation 3,415  to  3,660 
m. 

Wide  opening  in  the  heavily  tim- 
bered slope;  vertical  drop  390  m; 
length  823  m. 

Flat  bottom  of  the  Snake  River 
Valley;  under  severe  conditions  the 
avalanche  may  run  across  the  152- 
m-wide  valley  and  reach  the  high- 
way on  the  opposite  slope. 
The  starting  zone  is  controlled  in- 
tensively by  explosives  and  protec- 
tive skiing,  so  that  large  avalanches 
are  unlikely  unless  deposition  oc- 
curs from  a  prolonged  storm. 
By  a  process  described  in  detail  for 
the  Ironton  Park  avalanche  path, 
the  profile  for  the  Pallavicini  path  is 
computed  as  shown  in  figure  10.  Of 
particular  note  is  the  steep  80  m 
starting  zone,  in  which  the  avalanch- 
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ing  slab  initially  spans  the  upper  50 
m.  A  second  noteworthy  character- 
istic of  the  profile  is  the  adverse 
slope  beyond  the  creek  (cell  96)  of 
20%  average  grade. 

The  primary  consideration  in  evaluating  the  Palla- 
vicini  path  is  to  obtain  qualitative  information  on  the 
effect  of  the  adverse  grade  in  slowing  and  stopping 
avalanches.  Good  field  data  show  that,  in  this  path, 
only  deep  dry-snow  slab  avalanches  with  a  significant 
airborne  fraction  ever  reach  the  road.  Most  small 
midwinter  avalanches  stop  between  the  grade  change 
at  about  400  m  (cell  48)  and  the  creek  (cell  96),  (fig. 
10).  The  path  is  divided  into  10-m  increments  or 
cells,  and  the  avalanche  starting  zone  is  assumed  to 
be  50  m  long.  Table  3  summarizes  the  computer  runs 
made  using  the  Pallavicini  profile. 

An  avalanche  of  0.75  m  thickness  stops  quickly 
when  vicosity  and  friction  are  set  equal  to  0.55.  What 
might  be  considered  nominal  flow  for  an  avalanche 
0.75  m  in  thickness  occurs  when  viscosity  and  friction 
are  reduced  to  0.5.  In  this  case,  the  debris  stops  in 
cell  53,  which  is  a  bench  in  the  lower  part  of  the  track 
(fig.  10).  Flow  of  a  2-m  avalanche  with  v  and  f  at  0.5 
(table  3,  Run  8)  ends  40  m  beyond  the  creek,  but  does 
not  reach  the  road.  It  is  necessary  to  drop  v  and  f  to 
0.4  for  a  2-m  avalanche  to  reach  the  road  (Run  10). 
This  indicates  that  a  deep,  dry  avalanche  with  an 
airborne  phase  is  modeled  for  v  and  f  in  the  range  of 
0.4.  Since  many  avalanches  with  an  average  fracture 


face  height  of  2  m  on  the  Pallavicini  Path  are  known 
to  stop  in  the  vicinity  of  the  creek,  a  value  of  0.6 
seems  reasonable  for  v  and  f. 

The  20%  adverse  slope  beyond  the  creek  in  the 
Pallavicini  Path  is  an  effective  mechanism  for  stop- 
ping avalanches  that  flow  into  this  region,  and 
mounding  of  the  debris  is  noted  in  the  AVALNCH 
prediction  (table  3).  Reynolds  number,  Re,  which 
accounts  for  both  inertial  and  viscous  effects,  is 
computed  at  the  point  of  first  occurrence  of  zero 
slope  in  the  avalanche  path,  and  is  plotted  in  figure  1 1 
versus  runout  distance  for  two  of  the  paths  analyzed. 
For  Ironton  Park,  runout  is  across  a  frozen  lake  of 
zero  slope.  The  20%  adverse  slope  of  the  Pallavicini 
path  results  in  a  3.5  reduction  of  runout  distance  for 
Re  =  240,  and  up  to  a  5.7  reduction  at  Re  =  560. 
These  results  are  based  upon  a  nominal  flow  depth  of 
H  =  2  m.  Reynolds  number,  Re,  defined  as 

UH* 

Re  =   

V 

was  computed  using  U  =  the  velocity  of  the  ava- 
lanche leading  edge  at  the  instant  it  enters  the  runout 
zone,  H*  =  4H  and  v  is  the  kinematic  viscosity 
which  is  usually  near  0.5. 

So  far  we  have  concentrated  on  the  large  avalanche 
that  flows  onto  the  adverse  grade.  For  the  smaller 
avalanche  that  stops  in  the  track,  flow  height  is  a 
sensitive  parameter  that  controls  the  distance  of 
flow.  Comparing  Runs  1  and  2  in  table  3,  the 
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Table  3.  —  Travel  distance  for  the  Pallavicini  Avalanche  for  specified  flow  conditions 


Initial 

Run 

Viscosity 

Friction 

Thickness 

Travel 

Comments 

No. 

(mVs) 

Coeff. 

of  slab 

Distance 

(m) 

(m) 

0.55 

0.55 

0.75 

130 

Avalanche  stops  in  cell  17.  Maximum  debris  depth  is  0.7  m  and 

covers  80  m  of  slope. 

2 

0.50 

0.50 

0.75 

440 

Avalanche  stops  in  cell  53.  Maximum  debris  depth  is  0.6  m  and 

covers  130  m  of  slope. 

3 

0.5 

0.5 

1.00 

780 

Avalanche  stops  at  cell  90.  Maximum  depth,  1 .0  m;  extent  100  m. 

Initial  slab  40  m  extent. 

4 

0.6 

0.6 

1.00 

710 

Avalanche  stops  in  cell  82.  Maximum  debris  depth  is  0.8  m  and 

covers  170  m  of  slope. 

5 

0.5 

0.5 

1.50 

860 

Avalanche  stops  at  cell  97, 10  m  beyond  creek  maximum  depth  3.0 

m,  extent  40  m. 

6 

0.6 

0.6 

1.50 

850 

Avalanche  stops  at  creek  (cell  96).  Maximum  debris  depth  is  1 .5  m 

and  covers  60  m  of  slope. 

7 

0.6 

0.6 

2.00 

870 

Avalanche  stops  at  cell  98  just  beyond  creek.  Maximum  debris 

depth  7.4  m.  Avalanche  covers  40  m  of  slope. 

8 

0.5 

0.5 

2.00 

890 

Avalanche  stops  at  cell  100.  Maximum  depth  3.6  m,  extent  40  m. 

Adverse  slope  is  effective. 

9 

0.45 

0.45 

2.00 

910 

Avalanche  stops  at  cell  102.  Maximum  depth  9.0  m;  extent  40  m. 

10 

0.4 

0.4 

2.00 

Avalanche  leaves  grid  at  cell  105  with  velocity  of  5.0  m/s. 

sensitivity  of  travel  distance  to  the  values  of  v  and  f 
is  pronounced.  This  indicates  that  the  mechanisms  to 
stop  the  small  avalanche  on  a  positive  slope  are  the 
surface  friction  and  kinematic  viscosity  coupled  with 
local  variation  in  slope  profile.  Since  friction  and 
viscosity  are  not  well  defined  for  different  snow 
conditions,  it  is  unlikely  that  programming  for  the 
small  avalanche  will  be  systematically  correct  for  all 
avalanche  incidents. 

The  major  conclusions  we  draw  from  this  analysis 
of  the  Pallavicini  avalanche  path  are: 


1.  For  the  deep- slab,  high  air  content  avalanche 
flow,  values  of  v  and  f  of  the  order  0.4  yield 
results  consistent  with  observed  avalanche  be- 
havior on  this  slope. 

2 .  The  dominant  characteristic  of  adverse  slope  in 
halting  avalanches  is  noted,  which  substanti- 
ates long  standing  physical  evidence.  The  ef- 
fectiveness of  adverse  slopes  of  different  angles 
and  heights  to  stop  avalanches  is  considered  in 
greater  detail  in  the  section  headed  "  Sensitivity 
of  Runout  Distance  to  Various  Parameters." 


800 
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Hematite  Gulch  Avalanche  Path 


Location: 

Starting  Zone: 
Track: 


Runout  Zone: 
Avalanche: 


Profile: 


Highway  110  north  of  Howards- 
ville,  Colo. 

South  facing  bowl-shaped  basin. 
Vertical  drop  approximately  915 
m;  length,  approximately  670  m. 
hi  upper  reach,  flow  is  open  sheet 
flow,  which  changes  to  gully  flow 
over  55%  of  total  path. 
Broad  valley  of  bare  ground  and 
rocks  cut  by  the  Animas  River. 
A  wet  slab  avalanche  occurred  on 
May  7,  1975.  Fracture  thickness 
was  between  1  and  1.5  m  and 
tapered  off  at  the  flanks  to  about 
0.5  m.  Approximately  one-f0urth 
of  the  snow  in  the  starting  zone 
was  released.  The  flow  emerged 
from  the  gully  and  flowed  into  the 
Animas  River. 

Howardsville  Quadrangle,  Colo- 
rado (7.5').  A  centerhne  profile  of 
the  path  is  shown  in  figure  12. 
Three  bends  of  45,  30,  and  35 
degrees,  starting  with  the  upper- 
most, are  indicated  on  the  profile. 
Extent  of  the  gully  is  indicated 
also. 


The  slope-parallel  cell  dimension  is  taken  as  20  m 
in  this  case,  resulting  in  a  total  of  112  cells  to 
represent  the  entire  path.  An  embankment  occurs 
upslope  from  both  the  road  and  the  river  which  can 
not  be  dupHcated  accurately  in  AVALNCH,  since 
actual  modeling  would  require  depositing  snow  in  the 
depressions  formed  by  the  banks.^  Instead  the  road 
is  shown  one  cell  wider  than  actual,  with  these  cells 
level  (g^  =  0),  and  the  friction  set  at  f  =  0.9,  a  large 
value.  The  wet  snow  of  this  avalanche  is  assumed  to 
have  high  viscosity  {v  =  0.55).  Friction  is  also 
assumed  large,  at  f  =  0.55,  over  the  entire  path,  with 
local  increases  at  the  bends  (f  =  0.65),  and  in  the 
runout  zone  where  the  flow  crosses  beire,  rocky 
ground  (f  =  0.9).  Initial  slab  thickness  is  1.5  m  over  a 
downslope  span  of  40  m  (two  cells). 

The  leading  edge  of  the  avalanche  stopped  in  the 
Animas  River  (cell  110),  the  debris  extended  80  m 
upslope  from  the  river,  with  a  maximum  height  of  0.8 
m.  Maximum  velocity  of  the  flow  after  leaving  the 
steep  starting-zone  was  26  m/s  in  cell  15.  The  flow 
encountered  bare  ground  at  cell  89,  just  below  the 
lowest  bend  of  the  path.  The  velocity  of  the  leading 
edge  at  this  point  was  10  m/s. 


^Recent  research  with  the  computer  code  SMAC,  a  more 
elaborate  version  of  SOLASURF,  indicates  that  it  will  handle 
this  condition. 
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Figure  12.— Longitudinal  profile  of  Hematite  Gulch  avalanche  path 
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At  Hematite  Gulch  avalanches  are  known  to  cross 
the  road  and  river  and  flow  across  the  valley.  Increas- 
ing initial  slab  thickness  to  2  m  and  dropping  v  and  f 
values  of  0.5  and  0.4  approximately  models  the  worst 
case  avalanche  condition.  However,  specific  data  on 
the  worst  case  event  is  lacking  for  this  path,  so  no 
attempt  is  made  at  computer  simulation. 

Conclusions  drawn  from  this  analysis  are: 

1.  Viscosity  and  friction  coefficient  values  on  the 
order  of  0.55  appear  to  model  wet  spring  snow 
of  the  high  viscosity  type. 

2.  Friction  set  at  0.9  for  bare  ground  with  rocks 
appears  to  be  consistent  with  the  nominal  value 
of  0.55  for  the  remainder  of  the  profile  in 
stopping  the  avalanche  as  was  observed. 

3.  For  the  volume  of  snow  set  in  motion,  the  river 
and  road  represent  an  imposing  barrier  for 
stopping  the  flow. 

4.  Some  mechanism  should  be  incorporated  in 
AVALNCH  to  account  for  contracting  and 
spreading  of  the  flow  when  local  geometry 
dictates  three-dimensional  flow  conditions. 


Stanley  Avalanche  Path 


Location: 


Starting  Zone: 


Track: 


Runout  Zone: 


Avalanches: 


Front  range  Colorado;  south  of 
Berthoud  Pass;  southeast  slope  of 
the  east  shoulder  of  Stanley 
Mountain. 

No.  1.  A  bowl-shaped  depression 
on  the  slope  above  timberline;  ele- 
vation 3,540  to  3,780  m. 
No.  2.  Wind  blown  snow  slab  along 
a  treeline  extending  upward  from 
3,470  m  elevation. 
Vertical  drop  730  m;  length  1,700 
m.  U.S.  Highway  40  crosses  at  two 
elevations.  One  crossing  is  at 
about  1,300  m,  the  second  at 
about  1,700  m,  in  the  valley 
bottom. 

Lower  section  of  the  track  with  a 
gentle  slope  and  the  valley  of  Clear 
Creek. 

Stanley  frequently  blocks  the  up- 
per highway.  Most  avalanches  do 
not  reach  the  valley  bottom;  how- 
ever, a  few  of  the  larger  avalanches 
have  been  known  to  overrun  the 
lower  highway  and  come  to  rest  at 
the  foot  of  the  opposite  slope. 
Most  avalanches  releasing  at 
starting  zone  No.  2  stop  at  or 


above  the  upper  highway.  Most 
large  avalanches  released  in  start- 
ing zone  No.  1  stop  between  the 
two  highways. 
Profile:  Berthoud  Pass  Quadrangle,  Colo- 

rado (7.5').  A  centerline  profile  of 
the  eastern  most  of  the  three 
tracks  of  Stanley  is  shown  in  figure 
13.  The  slope  above  the  upper 
road  is  uniform  except  for  a  small 
rise  at  the  starting  zone. 


Major  consideration  is  modeling  the  upper  road  to 
represent  the  different  types  of  flow  adequately.  The 
upper  road  is  three  lanes  wide  with  a  turnout,  which  is 
equivalent  to  one  cell  of  20  m.  Two  additional  cells 
are  placed  at  zero  grade  to  represent  the  uphill 
highway  cut  and  the  snowplow  embankment  of  the 
downhill  side,  so  that  three  zero  grade  cells  (g^  =  0) 
constitute  the  road  intrusion.  These  three  cells,  plus 
one  additional  on  each  side  of  the  road,  are  specified 
with  large  friction  (f  =  0.9)  to  simulate  the  irregular 
snow  profile  in  the  vicinity  of  the  road.  For  the 
remainder  of  the  path,  typical  midwinter  hard  pack 
snow  coefficient  values  are  assumed  with  v  =f= 
0.55,  which  are  also  assumed  across  the  lower  road 
since  it  offers  little  resistance  to  flowing  snow. 

Results  of  several  "typical"  avalanche  runs  are 
listed  in  table  4.  We  see  that  1-m  deep  avalanches 
starting  in  release  zone  No.  2  are  trapped  at  the  road 
consistent  with  a  number  of  sightings.  A  1.5-m  deep 
avalanche  starting  in  release  zone  No.  1  barely  passes 
the  road,  then  regains  speed  and  finally  stops  be- 
tween the  roads.  Increasing  the  avalanche  thickness 
to  2  m  and  starting  in  release  zone  No.  1,  the 
avalanche  crosses  the  lower  road.  Tapered  avalanche 
slabs  2.5  and  3  m  thick  at  the  crown,  and  1  m  thick  at 
the  toe  100  m  downslope,  behave  similarly  to  the  2-m 
avalanche.  The  3-m  avalanche  flows  out  of  the  grid 
but  is  slowing  rapidly  on  the  adverse  grade  past  the 
creek.  Comparing  runs  Nos.  1  and  3,  and  Nos.  5  and 
6,  we  conclude  that  length  of  the  starting  slab  does 
not  change  flow  Hmits  significantly.  Flow  length  could 
be  important  when  the  leading  part  of  the  flow  fills  a 
recess  (road  cut,  creek  recess,  etc.)  and  the  flow  is 
long  enough  so  that  the  trailing  part  flows  over  the 
snow-filled  recess.  This  phenomenon  cannot  be 
modeled  exactly  by  AVALNCH  except  to  fill  the 
recess  by  changing  the  profile  of  the  path  prior  to 
making  the  computer  run. 

A  field  observation  on  March  4,1977,  showed  that 
an  avalanche  released  by  explosives  in  starting  zone 
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Figure  13.— Longitudinal  profile  of  Stanley  avalanche  path. 

Table  4.— Travel  distance  and  maximum  velocity  for  the  Stanley  Avalanche  for  specified  flow  conditions 


2000 


Initial 

Run 

Starting 

Viscosity 

Friction 

Thickness 

Maximum 

Cell 

Travel 

Cell 

No. 

Zone  No. 

(mVs) 

Coeff. 

of  slab 

Velocity 

No. 

Distance 

No. 

Comments 

(m) 

m/s 

(m) 

1 

1 

0.55 

0.55 

1.0 

9 

35 

1,300 

74 

Initial  40-m  slab  stops  in  upper 

road :  1 .3  m  max,  depth ;  extent 

1 00  m :  upper  road  f  =  0.55, 

2 

2 

0.60 

0.60 

1.0 

8 

34 

770 

73 

Initial  40-m  slab  stops  in  upper 

road;  1.1  m  max.  depth;  extent 

160  m;  upper  road  f  =  0.6. 

3 

2 

0.55 

0.55 

1.0 

10 

34 

800 

74 

Initial  100-m  slab  stops  in  upper 

road:  1 .8  max.  depth;  extent  280 

m:  upper  road  f  =  0.55. 

4 

1 

0.55 

0.55 

1.5 

22 

10 

1,600 

90 

Initial  60-m  slab  stops  between 

roads;  0.5  m  max.  depth;  extent 

320  m:  upper  road  f  =  0,9 

velocity  across  road  V  =  5  m/s, 

then  to  V  =13  m/s. 

5 

1 

0.55 

0.55 

2,0 

28 

10 

1,780 

99 

Initial  60-m  slab  stops  past 

second  road;  7,3  m  max,  depth: 

extent  20  m;  upper  road  f  =  0.9; 

velocity  across  road  V  =  23  m/s. 

6 

1 

0.55 

0.55 

2.0 

28 

74 

1,780 

99 

Initial  100-m  slab  stops  past 

tapered  to 

second  road;  8.6  m  max.  depth; 

1,0 

extent  20  m;  upper  road  f  =  0.55; 

velocity  across  road  V  =23m  s. 

7 

1 

0.55 

0.55 

2.0 

28 

74 

1 ,780 

99 

Initial  60-m  slab  stops  past 

second  road;  8.6  m  max.  depth; 

extent  20  m;  upper  road  f  =  0.55; 

velocity  across  road  V  =  26  m  s. 

8 

1 

0.55 

0.55 

2.5 

45 

75 

1,820 

101 

Initial  100-m  slab  stops  past 

tapered  to 

second  road;  9.0  m  max  depth; 

1.0 

extent  20  m;  upper  road  f  =  0.55, 

velocity  across  road  V  =  45  m/s. 

9 

1 

0.55 

0,55 

3.0 

60 

78 

Flows  out  of  grid 

Initial  100-m  slab  flows  out  of 

grid  at  V  =  37  m  '  s.  upper  road 
f  =0.55;  velocity  across  road 
V  =59m;s. 
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No.  1  reached  a  speed  of  33  m/s  above  the  upper 
road.  Most  debris  stopped  in  the  upper  road  (cell  70- 
73)  with  a  small  amount  extending  to  cell  No.  85. 
Average  fracture  face  thickness  was  1.51  m,  and 
average  slab  thickness  in  the  starting  zone  was  1.19 
m.  The  case  that  best  matches  this  occurrence  is  Run 
4  in  table  4.  The  greatest  discrepancy  is  in  the  two 
velocity  values,  which  is  probably  due  to  the  differ- 
ence between  the  slab  thickness  used  in  the  calcula- 
tions and  the  flow  height  actually  attained  in  the  field. 

The  upper  road  intrusion  modeled  for  Stanley  is 
apparently  sufficient  to  stop  the  smaller  avalanches, 
and  yet,  permit  the  larger  ones  to  flow  beyond.  User 
judgment  is  warranted  in  setting  up  these  intrusion 
models  for  each  case  that  is  encountered. 

Finally,  it  is  noted  that  on  the  near  uniform  slope  in 
release  zones  Nos.  1  and  2,  avalanche  velocity  for  the 
1-  and  1.5-m  avalanches  approximates  the  quoted 
condition  that  flow  velocity  reaches  80%  of  maximum 
speed  after  the  avalanche  has  covered  a  distance 
equal  to  25  times  the  initial  slab  thickness  (Voellmy 
1964).  An  alternate  statement  that  is  less  precise  is 
that  the  avalanche  accelerates  rapidly  at  the  start, 
then  flows  on  at  a  near  equilibrium  flow  rate  if  the 
slope  is  constant.  This  physical  condition  is  generally 
observed  in  the  AVALNCH  results. 

General  conclusions  derived  from  the  Stanley  ava- 
lanche study  include: 

1.  Heavy  midwinter  dry  snowpack  of  1-  to  2-m 
thick  slabs  may  be  modeled  with  v  =  f  =  0.55. 

2.  Slab  length  has  secondary  effect  on  avalanche 
flow,  under  restrictions  of  modeling  by 
AVALNCH. 

3.  Individual  road  cut  modeling  must  be  subject  to 
user  judgment.  The  specific  modeling  of  the 
upper  road  on  Stanley  as  outlined  herein  has 
yielded  physically  reasonable  results. 


Number  3  East  Avalanche  off  Max*  s  Mountain, 
Alyeska  Ski  Area,  Girdwood,  Alaska 


Location: 


Starting  Zone: 


Track 


Northwest  slope  of  Baumann 
Bump  (Max's  Mountain);  5  km  east 
of  Girdwood,  Alaska,  or  64  km 
southeast  of  Anchorage,  Alaska,  on 
the  north  side  of  Turnagain  Arm. 
One  of  a  series  of  narrow  rocky 
slots  at  an  elevation  of  about  760 
m. 

Open  slope  with  two  benches:  one 
at  about  400  m,  the  other  at  about 
230  m. 


Runout  Zone:      Deeply  incised  stream  valley. 

Avalanche:  Avalanches  often  start  in  cold,  dry 

snow  at  the  top  of  the  mountain 
and  encounter  damp  or  wet  snow 
on  the  lower  slopes.  Debris  stops  in 
the  stream  channel  with  only  air- 
borne particles  crossing  the  can- 
yon. 

Profile:  Seward  (D-6)  Quadrangle,  Alaska 

(15'  x  22.5').  A  centerline  profile  of 
the  No.  3  East,  Max's  Mountain 
path  is  shown  in  figure  14.  The 
several  rapid  changes  in  elevation 
on  the  profile  were  not  modified  in 
order  to  see  if  AVALNCH  could 
handle  them,  even  though  such 
rapid  changes  are  physically  un- 
likely because  of  snow  cover 
smoothing. 

The  avalanche  path  is  modeled  by  129  slope 
parallel  cells,  each  10  m  long.  A  frequent  condition  on 
this  path  is  dry  snow  down  to  an  elevation  of 
approximately  400  m  and  wet,  viscous  snow  below 
this  elevation.  All  avalanches  observed  in  the  past 
several  years  have  stopped  on  the  slope  or  in  the  gully 
at  cell  116  with  no  flow  onto  the  bench  of  the  gun 
tower. 

Three  runs  were  made  to  model  different  condi- 
tions known  to  occur  on  this  slope.  In  two  runs  the 
initial  slab  was  taken  as  1.5  m  thick  and  40  m  in 
extent.  In  the  third  run  slab  thickness  was  reduced  to 
1  m. 

Run  No.  1:  High  viscosity  (v=0.7)  is  assigned  the 
slab  to  simulate  wet  snow  in  the  start- 
ing zone.  Friction  is  set  at  f=0.7  to  cell 
80  and  ^0.9  to  cell  129,  to  model 
moderate  friction  on  the  upper  path 
and  large  friction  on  the  lower  path. 

Results  show  maximum  velocity  at- 
tained is  20  m/s  at  cell  12,  and  the 
avalanche  enters  cell  115  at  a  speed  of 
10  m/s  and  stops  in  cell  116.  The 
debris  is  scattered  over  120  m  to  a 
depth  of  0.94  m  in  cell  116. 
Run  No.  2:  Midwinter  deep  snowpack  viscosity  is 
assumed  with  v  =0.55  m^/s.  Friction 
is  set  atf=0.5  to  cell 72,  ^0.7  from  cell 
73  through  80,  and  ^0.9  for  the  re- 
mainder of  the  path.  This  is  intended  to 
simulate  dry  snow  running  on  a  dry 
upper  path  and  onto  wet  snow  on  the 
lower  path.  The  avalanche  attains  a 
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maximum  speed  of  27  m/s  at  cell  33, 
enters  cell  115  at  10  m/s,  and  stops  at 
cell  116.  Debris  covers  90  m,  to  a 
maximum  depth  of  1.4  m  in  cell  116. 
Run  No.  3:  The  same  viscosity  and  friction  condi- 
tions as  in  Run  No.  2  are  assumed.  The 
only  difference  is  that  the  thickness  of 
the  initially  released  slab  is  set  at  1  m. 
Under  these  conditions,  the  flow  does 
not  reach  the  creek,  but  stops  at  cell 
100.  As  noted  from  figure  14,  cell  100 
is  in  the  middle  of  a  flatter  part  of  the 
slope  that  starts  at  ceU  92. 

This  reduction  in  slope  coupled  with 
the  assigned  viscosity  and  friction  for 
this  region  are  apparently  sufficient  to 
halt  the  avalanche. 

Conclusions  drawn  from  this  example  are: 

1.  Friction  values  in  the  range  of  0.7  to  0.9  are 
large,  and  seldom  should  values  larger  than  this 
be  used  for  clear  slopes. 

2.  Program  AVALNCH  remains  stable  for  rapid 
variations  in  g^  when  the  slope  changes  sig- 
nificantly between  the  adjacent  cells.  Sign 
changes  for  g  in  the  case  of  adverse  slopes  are 
readily  handled  by  AVALNCK 
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No.  3  East  Max's  Mountain  path  is  a  good 
example  that  there  is  a  physical  basis  for 
treating  friction  as  a  variable  quantity  over  the 
path. 


Imogene  Avalanche  Path 


Location: 


Starting  Zone: 


Track: 

Runout  Zone: 
Avalanches: 


Profile: 


Second  unnamed  drainage  north  of 
middle  fork  of  Mineral  Creek;  top 
elevation  3,780  m;  east  facing 
slope. 

Broad  basin  with  large  rock  out- 
crops; lower  portion  has  small  rock 
outcrops  and  grass. 
Narrowing  V-shaped  depression, 
grass,  and  willows. 
Broad  fan,  grass  and  willows. 
Runout  is  into  the  Mineral  Creek 
valley.  An  avalanche  has  crossed 
the  highway  twice  since  1951.  Both 
runs  were  full-track.  Total  vertical 
drop  is  670  m;  length  of  path  is 
1,730  m. 

SUverton  Quadrangle,  Colorado 
(7.5').  This  path  is  characterized  by 
a  constricted  region  from  520  to 
880  m  distance  along  the  path  as 
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Figure  14.— Longitudinal  profile  of  Numbers  East  avalanche  path  (Max's  Mountain). 
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noted  in  figure  15.  After  flow 
crosses  the  creek  and  road,  the 
avalanche  encounters  adverse 
grade. 

The  path  is  modeled  by  91  cells  each  20  m  long. 
Initial  flow  height  is  set  at  1 .5  m,  and  the  extent  of  the 
starting  slab  along  the  fall  line  is  100  m.  Viscosity  is 
set  at  a  nominal  midwinter  value  of  0.5.  The  friction 
coefficient  is  varied  in  three  ways  in  the  constricted 
region  of  the  flow,  for  the  three  computer  runs  carried 
out. 


Run  No.  1: 


Run  No.  2: 


800r- 


Friction  coefficient  at  0.5  over  the 
entire  path.  Flow  velocity  is  maximum 
at  26  m/s  (cell  16)  and  decreases  until 
flow  stops  at  cell  85  midway  between 
the  creek  and  the  road.  Debris  height 
reaches  a  maximum  of  5.5  m  for  a 
distance  of  140  m.  The  leading  edge 
velocity  at  exit  from  the  constricted 
region  is  20  m/s. 

The  friction  coefficient  is  increased  to 
0.6  for  cells  in  the  constricted  region. 
Flow  stops,  as  in  run  No.  1,  at  cell  85. 
Debris  depth  is  5.5  m  maximum  over 
160-m  extent.  Although  stopping  dis- 


tance is  equal  for  the  two  cases,  the 
velocity  of  the  leading  edge  emerging 
from  the  constricted  region  in  this  case 
is  17  m/s,  a  decrease  of  159c  from  run 
No.  1. 

Run  No.  3:  The  friction  coefficient  is  decreased  to 
0.4  for  cells  in  the  constricted  region. 
Flow  again  stops  in  cell  85  with  maxi- 
mum debris  depth  of  5.6  m,  extending 
for  140  m.  Flow  velocity  coming  out  of 
the  constricted  region  is  about  24  m/s, 
an  increase  of  20%  over  that  of  run  No. 
1. 

Conclusions  drawn  from  the  study  of  this  avalanche 
path  are: 

1.  Varying  the  friction  coefficient  in  the  constric- 
ted zone  of  the  avalanche  track  had  no  effect  on 
total  runout  distance.  This  implies  that  travel 
distance  on  the  Imogene  path  is  dominantly 
controlled  by  the  geometry  of  the  runout  zone 
itself.  Generally,  the  dynamics  of  constrictions 
suggests  that  the  surface  friction  decreases  in 
importance  because  flow  depth  increases,  and 
viscosity  increases  in  importance  because  cir- 
culation increases.  However,  little  information 
exists  on  flow  of  snow,  so  that  we  have  no  basis 
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Figure  15.— Longitudinal  profile  of  Imogene  avalanche  path. 
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for  selecting  any  one  criterion  over  another.  Our 
previous  conclusions  that  v  and  f  have  approx- 
imately equal  influence  suggests  that,  in  the  case 
of  a  constriction,  no  adjustment  should  be  made. 

2.  Little  information  is  known  on  flow  properties 
of  snow,  and  further  research  is  warranted  in 
this  area. 

Numerical  Studies  with  Program  AVALNCH 

Equilibrium  Flow  Calculation 

To  establish  upper  limits  on  the  flow  velocity  that 
can  be  attained  as  a  function  of  slope  steepness, 
equilibrium  flow  calculations  were  carried  out  for  30-, 
40-,  and  50- degree  slopes.  Surface  friction  and 
viscosity  coefficients  were  assigned  values  from  0.4 
to  0.6,  which  is  representative  of  the  range  of  these 
parameters  in  avalanche  cases  considered  previous- 
ly. To  obtain  steady- state  flow,  continuity  in- flow  and 
out-flow  boundary  conditions  were  imposed  at  the 
upper  and  lower  boundaries,  respectively.  A  single 
cell  depth  of  flow  over  the  entire  domain  was  then 
specified  and  iteration  to  steady  flow  was  carried  out. 
Results  of  these  calculations  are  shown  in  figure  16, 
in  which  equilibrium  velocity  is  plotted  versus  flow 
height  of  the  core  material  with  slope  angle,  friction, 
and  viscosity  as  parameters.  In  actual  avalanche 
runout  cases,  the  flow  is  rotational  because  con- 
tinuous geometry  changes  cause  local  circulation. 
For  a  flow  to  achieve  the  velocities  shown  in  figure  16, 
it  is  necessary  for  the  geometry  to  be  constant  and 
continuous,  a  condition  that  is  unlikely  in  actual 
cases.  Thus,  the  velocities  shown  in  figure  16  can  be 
considered  upper  limits,  which  are  reported  only  for 
reference  when  interpreting  actual  flow  results. 

It  is  noted  that  in  AVALNCH  when  computing 
equihbrium  flow,  the  use  of  noninteger  heights  (e.g. 
1.5  m)  results  in  an  unstable  condition  that  is 
attributable  to  numerical  round- off  in  certain  in- 
structions in  the  code.  By  recoding,  the  problem  can 
be  averted;  however,  since  only  transient  conditions 
will  be  seen  in  all  avalanche  evaluations,  this  recoding 
is  not  attempted  in  the  unmodified  version  of 
AVALNCH. 

Variable  Width  Option  for  Program  AVALNCH 

The  width  of  many  avalanche  paths  vary.  The 
version  of  AVALNCH  used  previously  in  this  report 
does  not  account  for  variable  width.  The  interpreta- 
tion given  previously  in  the  results  section  should  be 
per-unit-width,  down  the  centerline  of  the  flow  as 
represented  by  the  longitudinal  profile  with  no  provi- 
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Figure  16.  — Equilibrium  velocity  as  a  function  of  flow  depth 
for  different  values  of  surface  friction,  f,  viscosity,  l',  and 
slope  angle,  (^. 

sions  made  for  spreading  or  narrowing  of  the  flow 
(two-dimensional  version  of  code).  For  avalanche 
paths  with  only  small  variations  in  the  flow  width,  the 
two-dimensional  formulation  is  adequate.  However, 
in  this  section  we  report  the  findings  of  a  pseudo 
simulation  of  three-dimensional  flow  in  which  the 
two-dimensional  flow  height  of  each  cell  is  adjusted 
each  CYCLE  of  iteration  to  account  for  local  changes 
in  flow  width.  This  approximation  is  not  a  complete 
three-  dimensional  formulation,  for  which  height  vari- 
ation along  a  contour  would  be  admitted.  However,  it 
will  account  for  large  variations  in  width  consistent 
with  the  accuracy  to  be  expected  based  upon  the 
nominal  slope  and  material  properties  we  are  able  to 
establish  for  flowing  snow. 

The  version  of  AVALNCH  including  the  variable 
width  option  is  listed  in  Appendix  D.  Arrows  are  used 
to  identify  changes  and  additions  from  the  nonvari- 
able  width  version  hsted  in  Appendix  B.  The  specific 
type  of  change  or  addition  can  be  determined  by 
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comparing  the  two  listings.  For  the  variable  width 
program,  additional  input  includes  a  width  factor, 
W(I),  for  each  cell  that  can  be  read  directly  from  a 
topographic  map.  The  width  factors  are  an  array  of 
numbers  that  follow  card  4  in  Appendix  A,  and  have 
the  same  format  as  card  4.  One  approach  to  setting 
values  for  is  to  select  a  value  of  1  for  cells  in  the 
starting  zone  of  the  avalanche  and  adjusting  this 
value  for  downslope  cells  as  the  width  changes.  In 
representing  width  change,  even  abrupt  changes,  the 
factors  should  vary  gradually,  since  snow  entrapment 
will  smooth  any  irregular  boundary. 

For  AVALNCH  to  accomodate  increases  in  height 
of  the  flow  due  to  a  constriction,  it  is  necessary  to 
include  more  than  one  row  of  cells  above  the  row  in 
which  the  avalanche  initially  flows. 

Appendix  D  provides  a  grid  five  cells  high  for  a 
flow  that  was  initially  only  one  cell  high.  This  re- 
quires a  total  of  seven  cells  high  (five  active  and  two 
boundary),  in  place  of  four  cells  of  the  non-variable 
height  code.  This  change  increases  the  storage 
allocation  of  the  code,  and  computer  running  time. 
For  the  Imogene  path,  described  below,  computer 
time  is  more  than  doubled  to  account  for  variable 
width. 


The  Lnogene  avalanche  path  is  used  as  an  example 
of  a  path  with  variable  width.  The  outhne  of  this  path 
is  depicted  by  the  solid  lines  in  figures  17  and  18.  The 
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Figure  17.— Topography  map  of  Imogene  avalanche  path, 
from  which  geometric  path  data  was  taken. 


Figure  18.— Photograph  of  Imogene  and  adjacent  ava- 
lanche paths  (Miller,  Armstrong  and  Armstrong  1976). 


shaded  region  in  figure  17  is  the  assumed  extent  of 
the  starting  zone.  The  dashed  lines  are  the  assumed 
boundaries  of  the  flow  after  release.  Under  these 
assumptions,  the  width  factor  read  directly  from 
figure  17  is  plotted  in  figure  19,  from  which  indiv- 
idual cell  values  are  interpolated. 

Results  of  the  variable  flow  calculations  show  that 
the  flowing  snow  emerges  from  the  constricted  region 
of  the  path  at  a  speed  of  approximately  23  m/s  and 
comes  to  rest  in  cell  85.  Thus,  Run  No.  3  of  the 
Imogene  path  is  roughly  duplicated,  indicating  that 
by  decreasing  the  friction  coefficient  in  the  constrict- 
ed area,  the  overall  effect  of  flow  constriction  is 
modeled.  We  note  that  the  variable  width  flow  stops 
in  cell  85,  the  same  as  reported  in  the  section  on  the 
Imogene  path  for  constant  width  flows.  Cell  85  is  the 
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first  cell  of  adverse  grade  (20%)  after  the  flow  crosses 
the  creek  bed.  This  amount  of  adverse  grade  is 
significant  in  stopping  snow  flow  over  a  broad  range 
of  entering  velocities.  Thus,  a  meaningful  number  for 
comparison  of  the  different  runs  is  the  velocity  of  the 
leading  edge  of  the  flow  entering  cell  85.  The  leading 
edge  velocity  and  the  final  maximum  debris  height 
are  listed  in  table  5  for  the  variable  width  runs  of  this 
section  and  the  three  constant-width  runs  of  the 
Imogene  path  section.  The  reported  velocity  is  the 
average  of  the  velocities  computed  for  cells  80 
through  84,  which,  by  averaging,  reduces  total  com- 
putational variations. 

We  conclude  from  these  results  that  the  various 
assumptions  made  in  modeling  the  constricted  zone, 
bring  about  an  overall  effect  on  terminal  flow  that  is 
small  compared  to  the  influence  of  terminal  zone 
geometry.  Significant  variation  in  snow  debris  height 
is  noted  from  table  5.  For  the  constant  flow- width 
cases,  two  active  cells  are  used  to  represent  vertical 
motion.  By  comparing  the  two- cell- high  runs  with  the 
variable-width  runs  having  three,  four  and  five  cells 
vertically,  we  note  apparent  numerical  instability  in 
height  prediction  when  flow  enters  the  upper 
boundary  cell  layer.  The  same  effect  occurs  also  in 
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the  variable-width  runs  in  the  constricted  zone  where 
flow  height  increases  by  a  factor  of  2.5.  With  this 
increase  in  height,  flow  is  in  the  cell  next  to  the  upper 
boundary  cell  for  the  three-cell-height  case,  and  still 
one  cell  or  more  below  the  boundary  cell  in  the  four 
and  five  cell  cases.  The  subsequent  differences 
arising  from  interpolating  with  a  boundary  cell  in  one 
case,  and  with  an  active  cell  in  the  other,  is  the  source 
for  the  differences  in  terminal  speeds  between  the 
three-  and  four-cell  cases.  For  consistent  results  one 
cell  height  should  be  maintained  between  the  cell  of 
maximum  flow  height  and  the  upper  boundary  cell  at 
all  points  along  the  path.  It  appears  that  if  accurate 
prediction  of  debris  height  is  not  required,  all  models 
predict  terminal  velocities  within  the  limits  of  ac- 
curacy encountered  on  selecting  values  for  v  and  f 
(to  represent  material  and  surface  conditions). 

Based  upon  these  results,  we  conclude  that  we 
have  a  choice  of  ways  to  account  for  moderate 
constriction  (or  expansion)  of  factors  up  to  2  or  3. 
One  approach  is  to  use  the  variable-width  code  at  a 
cost  more  than  double  that  for  the  constant-width 
code.  Alternatively,  we  can  use  the  constant-width 
code  and  reduce  friction  in  the  constricted  region. 
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Figure  19.— Longitudinal  profile  of  Imogene  avalanche  path  and  width  fraction. 
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Table  5.— Summary  of  leading  edge  terminal  velocity  and 
maximum  debris  height  for  the  Imogene  Avalanche  path. 


Run 
No. 


Description  of 
numerical  simulation 
in  computer  run 


Average  flow 
velocity  into 
cell  85  (m/s) 


Maximum 

height  of 
snow  debris 
(m) 


1  Constant  width  flow. 
Friction  f  =  0.5  over  path, 

and  f  =0.4  in  constriction  10.98  5.6 

region.  Two-cell-height 
model  of  flow  regime. 

2  Constant  width  flow. 
Friction  F  =  0.5  over  entire 

path.  Two-cell-height  10.74  5.5 

model  of  flow  regime. 

3  Constant  width  flow. 
Friction  f  =  0.5  over  path , 

and  f  =  0.6  in  constriction  10.68  5.5 

region.  Two-cell-height 
model  of  flow  regime. 

4  Variable  width  flow. 
Friction  f  =  0.5  over  entire 

path.  Three-cell-height  11.9  2.9 

model  of  flow  regime. 

5  Variable  width  flow. 
Friction  f  =  0.5  over  entire 

path.  Four-cell-height  12.9  3.1 

model  of  flow  regime. 

6  Variable  width  flow. 
Friction  f  =  0.5  over  entire 

path.  Five-cell-height  12.9  3.2 

model  of  flow  regime. 


This  is  physically  incorrect  because  it  does  not 
account  for  proper  mass  distribution  but  it  gives 
acceptable  results.  If  prediction  of  the  distribution  of 
terminal  debris  in  cases  of  pile-up  on  adverse  grade  is 
required,  then  cost  differences  between  the  two 
methods  is  removed.  This  occurs  because  a  height  of 
more  than  two  active  cells  is  then  needed  in  the 
constant-width  code  also,  which  is  the  factor  that 
increases  cost. 

A  possible  appHcation  of  the  variable-width  code  is 
in  modeling  avalanche  flow  in  which  flow  from  a 
tributary  comes  into  that  of  the  main  stream.  In  this 
case,  the  tributary  inflow  is  estimated  in  terms  of  the 
fractional  increase  it  introduces  over  mainstream 
flow  height.  Denoting  this  factor  by  r)  ■  and  main 
stream  flow  width  fraction  by  Wj  for  cell  j  on  the  path, 
then  the  adjustment  of  flow  height,  accounting  for 
both  width  variation  and  tributary  inflow,  is 

W;    .  . 


'J, 


W; 


In  actual  application  of  this  equation,  account  must 
also  be  taken  of  the  time  sequencing  between  main 


and  tributary  flows.  This  feature  has  not  been  pro- 
grammed into  the  variable-width  code,  but  it  war- 
rants consideration  in  order  to  evaluate  the  worst- 
case  occurrence  for  some  types  of  avalanche  paths. 


Sensitivity  of  Runout  Distance  to  Various 
Parameters 

Results  are  reported  from  numerical  experimenta- 
tion carried  out  with  the  Ironton  Park  avalanche  path. 
Runout  of  this  path  is  onto  a  flat  frozen  lake  surface. 
Because  of  the  simple  geometry  of  this  path,  altera- 
tion of  the  path  was  introduced  into  the  computer 
model  to  assess  the  sensitivity  of  the  flow  to  different 
parameters.  Among  parameters  considered  were 
average  thickness  of  the  initially  released  slab,  fric- 
tion and  viscosity  of  the  flow,  adverse  grade  in  the 
runout  zone,  and  spreading  angle  of  the  flow  in  the 
runout  zone.  The  thickness,  friction,  and  viscosity 
parameters  have  already  been  evaluated  and  are 
reported  in  table  2.  The  adverse  grade  and  terminal 
flow  speading  angle  parameters  remain  to  be  as- 
sessed. These  were  introduced  separately  into  the 
Ironton  Park  runout  by  artificially  representing  their 
effect.  For  example,  although  an  actual  avalanche 
would  not  spread  appreciably  on  the  flat  surface, 
numerical  spreading  was  introduced  in  the  computer 
model.  Thus,  the  effect  of  spreading  on  runout 
distance  was  appraised  without  introducing  the  at- 
tendant slope  rolloff,  which  must  actually  produce 
the  phenomenon.  The  effect  of  adverse  grade  was 
likewise  modeled  by  simple  specification  of  gravity 
components  in  place  of  those  corresponding  to  the 
lake  surface. 

In  the  modeling  of  flow  spreading,  angle  0  was 
introduced  starting  at  cell  73  (fig.  8),  which  is  the  first 
cell  on  the  surface  of  the  lake.  Angle  0/2  is  then  the 
angle  between  the  centerline  of  the  flow  and  one  of 
the  flow  boundaries.  Results  are  presented  in  figure 
20  of  spreading  angle,  0,  versus  runout  distance  for 
V  =  f  =  0.5  and  an  initial  release  slab  thickness  of  2 
m.  Spreading  begins  at  615  m  distance  (cell  73).  An 
angle  change  from  0  to  45  degrees,  holding  all  other 
parameters  constant,  decreased  travel  distance  by 
21%. 

Results  of  various  adverse  slopes  in  the  runout 
zone  of  the  Ironton  Park  path  are  shown  in  figure  21. 
As  expected,  travel  distance  continues  to  shorten  as 
the  adverse  grade  steepens,  finally  becoming  asymp- 
totic at  615  m  for  an  infinitely  steep  slope. 

Using  the  data  in  figures  20  and  21  and  table  2,  a 
measure  of  the  sensitivity  of  runout  distance  to  the 
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various  parameters  can  be  formulated  for  the  Ironton 
Park  avalanche  path.  For  example,  by  selecting  a 
nominal  set  of  values  for  the  different  parameters,  we 
can  then  vary  a  specific  parameter  by  a  given  amount 
and  note  the  change  in  runout  distance.  Spreading 
angle  was  nominally  set  at  0  =  10  degrees,  and  20% 
and  407c  increases  in  angle  were  introduced  to 
evaluate  change  in  runout  distance.  Adverse  slope 
was  nominally  set  at  10  degrees,  and  20%  and  40% 
increases  were  evaluated.  The  same  procedure  was 
followed  using  nominal  values  of  viscosity  and  fric- 
tion at  0.5,  and  slab  thickness  at  2  m.  Results  of  these 
computations  are  shown  in  figure  22.  It  should  be 
noted  that  most  of  the  variations  of  parameters  with 
runout  distance  are  nonlinear  (figs.  9,  20,  and  21), 
and  the  selection  of  a  nominal  value  was  arbitrary.  An 
attempt  was  made  to  select  values  that  did  not  lie  in 
either  a  flat  or  a  steep  portion  of  any  curve. 

The  intent,  then,  was  to  establish  order-of-mag- 
nitude  estimates  of  the  relative  importance  of  the 
parameters  for  a  simple  avalanche  path.  The  extent 
to  which  this  goal  was  achieved  is  indicated  in  figure 
22.  Runout  distance  is  most  sensitive  to  the  thickness 
of  slab  with  a  mean  effect  that  a  30%  variation  in  slab 
thickness  results  in  a  70%  change  in  runout  distance. 
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Figure  21.— Avalanche  runout  distance  as  a  function  of 
adverse  slope  angle  for  the  Ironton  Park  avalanche  path. 
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Figure  20.— Avalanche  runout  distance  and  maximum  ter- 
minal debris  height  as  a  function  of  flow  spreading  angle 
for  the  Ironton  Park  avalanche  path. 


Figure  22.— Change  in  avalanche  runout  distance  as  a 
function  of  change  in  various  flow  parameters  for  the 
Ironton  Park  avalanche  path. 


23 


This  is  followed  by  viscosity  or  friction  in  which  a  30% 
change  in  parameter  results  in  a  45%  change  in 
runout  distance.  Adverse  slope  change  of  30%  cor- 
responds to  25%  change  in  runout  distance,  and  30% 
change  in  spreading  angle  amounts  to  only  about  a 
2%  change  in  runout  distance.  One  parameter  not 
reported  in  figure  22  is  the  slope-parallel  length  of 
snow  slab  released,  but  it  is  determined  to  have 
negligible  effect  on  runout  distance  (see  the  section 
on  Stanley  path.) 

The  primary  conclusion  we  draw  from  these  results 
is  that  the  parameters  of  adverse  grade  and  spread- 
ing angle,  which  can  most  readily  be  evaluated  from 
topography  maps,  have  the  least  effect  on  runout 
distance.  The  parameters  of  slab  thickness  and 
friction  and  viscosity  of  the  flow,  have  the  greatest 
effect  on  runout  distance.  These  are  much  more 
difficult  to  evaluate,  and  certainly  deserve  special 
attention  in  future  avalanche  research. 


Leading  Edge  Snovv^  Entrainment  Modeling 

Many  avalanches  entrain  additional  snow  as  they 
move  down  the  slope.  Generally,  the  entrained  snow 
is  loose  and  lightweight — the  result  of  recent  snowfall. 
Thus,  we  envision  two  possible  influences  on  flow. 
First,  the  volume  of  the  moving  mass  increases;  since 
more  mass  is  being  accelerated,  flow  velocity  at  any 
point  on  the  path  should  be  lower  than  what  it  would  be 
in  the  absence  of  entrainment.  Second,  the  loose  snow 
at  the  leading  edge  may  reduce  friction  by  air 
entrainment  at  the  flow  boundary  interface.  This 
would  speed  up  the  flow.  To  model  the  case  of 
reduced  friction  in  AVALNCH,  we  simply  decrease 
the  value  of  the  friction  coefficient,  f.  We  noted 
previously  that  reducing  friction  to  f  =  0.4  on  the 
Pallavicini  path  results  in  snow  runout  that  duplicat- 
ed avalanches  with  a  significant  airborne  component. 

To  model  the  effect  of  an  increase  in  mass  due  to 
entrainment  requires  modifications  to  AVALNCH.  As 
currently  written,  program  AVALNCH  converts  ele- 
vation change  of  each  finite  difference  cell  along  the 
path  to  components  of  gravity,  g^  and  g^.  Hence, 
gravity  components  along  the  entire  path  are  speci- 
fied initially.  Thus,  if  we  attempt  to  simply  represent 
the  snow  mass  at  the  leading  edge  by  specifying  a 
non-zero  height  of  snow  in  all  cells,  then  all  snow  on 
the  sloping  path  will  accelerate  since  g  and  g^  are 
specified.  This  is  not  the  desired  result.  We  want  the 
avalanche  to  encounter  stationary  snow,  which  is 
then  brought  up  to  the  speed  of  the  advancing 
volume.  To  model  this,  consistent  with  the  finite 


difference  representation  used,  we  simply  introduce 
one  cell  of  stationary  snow  at  the  leading  edge  of  the 
avalanche  each  time  the  flow  enters  a  new  cell.  The 
momentum  balancing  the  ITER  phase  of  the  numer- 
ical methodology  in  AVALNCH  accelerates  this  new 
snow  and  assimilates  it  into  the  mass  of  flowing  snow, 
with  a  corresponding  increase  in  volume. 

The  following  changes  are  necessary  to  introduce 
constant- depth  leading- edge  snow  entrainment  over 
the  entire  length  of  an  avalanche  path: 

1.  Specify  thickness  of  entrained  snow  by  para- 
meter HNT,  and  insert  instruction  of  the  form 
HNT=  0.05  (units  of  meters)  following  instruc- 
tion 55  in  Appendix  B. 

2 .  Insert  the  following  statements  after  instruction 
264  in  Appendix  B: 

HIT  =  H(LDEG) 

H(LDEG)  =  HIT  +  HNT 

U(LDEG,  2)  =  HIT  *  U(LDEG,  2)/H(LDEG) 

U(LDEG,  1)  =  U(LDEG,  2)  *  (1.0—2.0  * 

FRC(LDEG)) 

V(LDEG,  2)  =  HIT  *  V(LDEG,  2)/H(LDEG) 

FV0L  =  FV0L  +  HNT  *  DELX 

FV0L1  =  FV0L1  +  HNT  *  DELX 

If  a  variable  entrainment  thickness  is  to  be  modeled, 
then  a  one-dimensional  array  must  be  specified  and 
the  thickness  values  read  into  AVALNCH  as  part  of 
the  input  data. 

As  a  simple  example,  the  Ironton  Park  avalanche 
path  described  previously  is  used  to  evaluate  the 
effect  of  mass  entrainment.  Additional  mass,  ex- 
pressed as  a  fraction  of  the  maximum  avalanche  slab 
thickness  of  2  m,  is  assumed  over  the  entire  path  (fig. 
8).  Based  upon  snow  entrainment  thickness  of  0.02  m 
(1%  of  slab  thickness),  0.05  and  0.1  m,  the  nominal 
leading  edge  velocity  versus  cell  number  along  the 
path,  is  shown  in  figure  23  for  each  case.  It  is  seen  that 
slowdown  occurs  as  entrainment  thickness  is  in- 
creased. However,  as  the  flow  comes  off  the  sloping 
part  of  the  path,  the  total  momentum  of  the  moving 
mass  must  increase  with  an  increase  in  entrainment 
thickness,  because  more  snow  is  accelerated  by 
gravity  than  when  entrainment  is  absent.  Thus,  if 
entrainment  were  eliminated  in  the  runout  zone,  the 
enlarged  avalanche  mass  (due  to  entrainment  on  the 
slope)  entering  the  runout  zone  should  flow  farther 
than  in  the  case  of  no  entrainment.  The  imposed 
continued  entrainment  in  the  runout  zone  acts  to 
resist  the  continued  motion  of  the  mass.  For  Ironton 
Park  the  two  effects  apparently  cancel,  and  runout 
distances  of  the  different  cases  do  not  vary  signifi- 
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■  No  snow  entrainment. 

60  -   Depth  of  entrained  snow  equal  to  1  %  of 

initial  slab  thickness. 

 Depth  of  entrained  snow 

equal  to  2.5%  of  initial  slab  thickness. 

 Depth  of  entrained  snow  equal  to  5%  of 

50  "  initial  slab  thickness. 


E 


CELL  NUMBER 

Figure  23.— Leading  edge  velocity  variation  with  cell  number  along  the  Ironton  Park 
avalanche  path  for  different  leading  edge  mass  entrainment  heights. 


cantly.  It  is  expected,  however,  that  for  paths  with 
short  runout  zones,  such  as  the  Pallavacini,  ava- 
lanches with  entrained  snow  would  run  further  than 
those  without. 

The  results  of  figure  23  are  based  upon  a  two-cell 
vertical  grid,  for  which  flow  heights  of  the  snow 
entrainment  cases  exceeded  values  that  correspond 
to  accurate  prediction  of  snow  distribution.  Thus,  it 
is  not  possible  from  the  data  obtained  to  compute 
momentum  values  for  the  different  cases.  However, 
based  upon  the  findings  in  the  section  headed 
"Variable  Width  Option  for  Program  AVALNCH," 
we  can  expect  the  predicted  velocities  and  travel 
distances  to  be  realistic  and  slow  down  in  the  runout 
zone  to  be  more  uniform  and  gradual  than  for  the  z  ero 
entrainment  case.  A  more  accurate  five- cell  vertical 
model  was  run  for  the  5%  entrainment  case  and 
showed  that  momentum  of  the  entrainment  case  does 
not  exceed  that  of  the  zero  entrainment  case.  It  also 
showed  runout  to  be  one  or  two  cells  farther  than  is 
reported  for  the  two- cell  vertical  case.  All  of  the 
entrainment  cases  could  be  rerun  using  a  more 
refined  vertical  scale  (with  increase  in  computer 
expense);  however,  the  results  presented  demon- 
strate that  there  is  no  fundamental  difficulty  in 


representing  the  entrainment  problem.  More  ration- 
al evaluation  can  be  carried  out  when  more  actual 
avalanche  events  are  available  for  reference. 

Discussion  and  Summary 

A  computer  code  applicable  for  predicting  snow 
avalanche  runout  distance  has  been  developed  based 
upon  two-dimensional,  viscous  fluid  dynamic  equa- 
tions. Basic  parameters  of  surface  friction,  fluid 
viscosity,  and  local  terrain  variations  are  accounted 
for  separately  by  empirical  evaluations  based  on 
known  avalanche  occurrences.  The  primary  objec- 
tive in  developing  the  code  was  to  model  gross 
avalanche  response  and  limiting-case  events.  Only 
gross  response  was  modeled  because  of  lack  of 
knowledge  about  the  properties  of  flowing  snow.  In 
particular,  specification  of  internal  viscosity  and 
surfac.e  friction  as  functions  of  distance  along  the 
avalanche  path  are  related  in  some,  as  yet  unknown, 
manner  to  particle  size  in  the  flow,  surface  condi- 
tions, air  entrainment,  density  distribution,  and 
other  factors  that  make  each  avalanche  distinct. 

In  modeling  gross  response,  local  variations  or 
gradual  changes  in  parameter  values  with  travel 
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along  the  path  are  assumed  to  be  averaged  out, 
provided  the  dominant  characteristics  of  the  flow  are 
represented.  Empirical  evaluations  are  based  on  spe- 
cific avalanche  events  in  order  to  successively  de- 
velop technology  in  modeling.  A  simple  avalanche 
slope  (see  section  on  Ironton  Park  Avalanche  Path) 
was  used  to  establish  nominal,  total-path  values  for 
the  coefficients  of  surface  friction  and  internal  vis- 
cosity based  upon  known  early  winter  avalanches.  It 
was  established  that,  with  the  units  of  measure  used, 
friction  and  viscosity  coefficients  of  approximately 
equal  value  and  of  the  order  f  =  u  =  0.5  model  the  flow. 
In  addition,  some  technique  is  needed  to  duplicate  the 
slowdown  dynamics  of  avalanches.  An  exponential 
dependence  of  surface  friction  on  avalanche  velocity 
at  the  time  of  terminal  flow  was  selected  based  on 
numerous  observations,  although  rare  exceptions 
are  noted  where  avalanches  have  flowed  inordinate 
distance  for  reasons  unknown. 

The  second  avalanche  path  studied  is  the  Palla- 
vicini,  where  a  terminal  adverse  grade  is  effective  in 
stopping  the  larger  avalanches.  On  this  path,  it  is 
known  that  only  avalanches  having  a  significant 
airborne  fraction  reach  the  road  (fig.  10).  To  match 
this  occurrence,  friction  and  viscosity  were  reduced 
to  V  =  f  =  0.4.  It  should  be  emphasizedthat  only  the 
more  consolidated  central  core  of  the  avalanche,  not 
the  turbulent  airborne  dust  cloud,  is  modeled.  Data 
are  scarce  concerning  airborne  avalanches  or  the 
airborne  component  of  mixed  motion  avalanches. 
There  is  even  a  lack  of  information  on  the  distribution 
of  terminal  debris  from  such  avalanches  which  would 
help  define  the  amount  and  distribution  of  snow  in 
the  snow  dust  cloud. 

The  third  and  fourth  examples  reported  are  the 
Hematite  Gulch  and  Stanley  avalanche  paths,  which 
are  considered  because  of  their  local  terrain  anoma- 
lies. Hematite  Gulch  has  three  bends  in  the  path,  and 
the  avalanche  modeled  ran  out  onto  rocky,  bare 
ground  (fig.  12).  Gross  modeling  of  this  spring,  wet 
snow  avalanche  was  achieved  by  letting  f=  v  =0.55, 
for  most  of  the  path,  with  values  of  f  =  v  =  0.65  at 
the  bends.  For  runout  onto  the  rocky  ground,  friction 
was  set  at  f  =  0.9  and  viscosity  left  at  v  =  0.55.  For 
these  values  the  model  showed  the  avalanche  flowing 
to  the  recorded  terminal  debris  zone.  Lacking  suffi- 
cient field  data,  further  evaluation  of  local  terrain 
effects,  as  depicted  by  the  Hematite  Gulch  path,  is 
not  warranted.  Flow  velocities  along  the  path  as  the 
avalanche  advances  are  additional  basic  data  needed 
to  see  how  well  the  model  duplicates  the  actual  event. 
This  type  of  data  would  result  in  more  accurate 
specifications  of  f  and  v  . 


The  Stanley  path  (fig.  13)  is  a  rather  uniform  slope 
interrupted  by  two  highway  cuts.  Small  to  moderate 
size  avalanches  are  known  to  be  trapped  by  the  upper 
road  cut,  and  large  avalanches  to  flow  to  the  lower 
road  cut.  Predicted  response  was  obtained  using 
midwinter  snowpack  conditions  v  =  f  =  0.55  over 
the  entire  path.  At  the  road  cuts,  friction  was  in- 
creased to  f  =  0.9  to  represent,  in  some  sense,  the 
resistance  of  the  road  bed  and  snow  embankments  to 
flow. 

The  fifth  avalanche  path  analyzed  is  at  Alyeska, 
Alaska,  very  close  to  the  ocean  (see  section  on 
Number  3  East  Avalanche).  The  typical  effect  of  the 
maritime  climate  is  a  wet,  viscous  snowpack  on  the 
lower  part  of  the  slope  and  drier  snow  farther  up  the 
mountain.  Thus,  the  initially  released  slab  is  often 
relatively  dry  {  v  =  0.55)  and  flows  over  a  similar 
snowpack  (f  =  0.5).  Farther  along  the  path,  friction  is 
increased  to  f  =  0.7  to  0.9  to  account  for  the  viscous 
effects  of  the  wet  snow  in  the  lower  part  of  the  path. 
Based  upon  this  modeling,  a  slab  1  m  deep  and  40  m 
long  in  the  slope-parallel  direction  stopped  short  of 
the  gully  at  the  foot  of  the  slope  (fig.  14),  whereas  a 
similar  slab  1.5  m  deep  flowed  into  the  gully.  Data 
from  the  area  are  rather  general,  however,  and  about 
all  we  can  say  is  that  v  and  f  values  are  consistent 
with  those  found  satisfactory  for  other  evaluations. 

The  final  example  reported  is  the  Imogene  path 
which  depicts  the  classic  "dog  bone"  or  contraction- 
expansion  flow  common  to  many  paths.  The  first 
attempt  to  model  this  variable  width  characteristic 
involved  adjustment  of  the  friction  coefficient,  f. 
However,  evidence  is  lacking  as  to  whether  f  should 
be  increased  or  decreased  particularly  for  a  viscous 
material  like  snow.  Later  a  pseudo-variable-width 
option  added  to  program  AVALNCH  verified  that  f 
should  be  increased.  This  result,  the  outcome  of  a 
numerical  analysis,  warrants  further  verification 
from  field  measurements. 

Having  established  a  basis  for  selection  of  values 
for  V  and  f  to  model  "typical"  cases  of  avalanche  flow, 
we  next  studied  particular  aspects  of  avalanche  flow 
to  which  our  numerical  model  is  applicable.  This  is 
not  to  suggest  that  specification  of  v  and  f  is  automatic 
for  different  avalanche  cases.  Each  avalanche  evalua- 
tion must  be  considered  seperately,  but  as  experi- 
ence is  gained  in  specifying  V  and  f  and  assessing  the 
results  of  analysis,  credibility  of  the  methodology 
should  improve.  One  factor  that  aids  in  establishing 
credibility  is  the  small  range  of  v  and  f  needed  to 
model  a  wide  variety  of  avalanches  from  the  large 
airborne  avalanche  to  the  wet,  viscous,  spring  ava- 
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lanche.  This  advantage  over  other  previously  devel- 
oped analysis  schemes  is  the  direct  result  of  separat- 
ing flow  properties  from  local  geometry  effects.  An 
additional  feature  of  program  AVALNCH  is  that  the 
basic  numerical  approach  of  subdividing  the  flow 
domain  into  an  array  of  cells  permits  the  user  to 
incorporate  additional  information  later. 

Since  open  channel  equilibrium  flow  has  been  used 
in  most  dynamic  evaluations  of  avalanches,  one 
numerical  study  using  program  AVALNCH  was  to 
establish  equilibrium  flow  velocity  as  functions  of  v, 
f,  flow  height,  and  slope  angle.  It  is  not  anticipated 
that  this  data  will  be  extremely  useful  in  the  future; 
however,  it  is  included  for  completeness. 

As  mentioned  above,  a  variable  width  option  of 
program  AVALNCH  is  reported  that  permits  evalua- 
tion of  snow  flow  as  the  width  of  the  flow  changes. 
This  phenomenon  occurs  when  snow  flows  into  a 
constricting  region  such  as  a  gully,  or  into  an  expand- 
ing region  such  as  an  alluvial  fan.  The  three-dimen- 
sional effect  is  represented  by  maintaining  continuity 
of  flow,  not  by  a  more  accurate  accounting  which 
would  require  the  introduction  of  a  three-dimension- 
al grid.  This  is  another  example  where  a  numerical 
capability  was  developed,  but  cannot  be  verified  until 
better  field  data,  including  flow  velocities,  are  avail- 
able. The  variable  width  runs  on  the  Imogene  path 
establish  the  friction,  f,  should  be  increased  in 
constricting  flow  when  the  constant  width  version  of 
the  code  is  used. 

Perhaps  the  most  significant  numerical  study  report- 
ed is  the  sensitivity  or  error  analysis  of  the  basic 
pararmeters  in  avalanche  flow.  These  parameters  are 
the  friction  and  kinematic  viscosity,  slab  thickness, 
adverse  grade,  and  flow  spreading  angle.  By  estabUsh- 
ing  nominal  values  for  these  parameters,  then  allowing 
one  parameter  to  vary,  the  change  in  avalanche  runout 
distance  can  be  computed.  For  nominal  30%  changes 
in  each  of  these  parameters,  imposed  separately, 
runout  distance  is  found  to  change  by  70%  for  slab 
thickness,  45%  for  friction  or  viscosity,  25%  for  adverse 
grade,  and  2%  for  spreading  angle.  Thus,  a  scale  of  the 
relative  importance  of  different  types  of  field  measure- 
ments and  observations  needed  for  a  better  definition 
of  snow  avalanche  dynamics  is  set. 

The  final  study  involves  evaluation  of  the  effect  on 
avalanche  dynamics  when  loose  snow  is  entrained  into 
the  flow  at  the  leading  edge  of  the  avalanche.  Leading 
edge  entrainment  can  result  in  change  of  the  surface 
friction,  as  well  as  gradual  increase  of  the  volume  of 


snow  in  motion.  In  the  study  reported,  only  the  effect 
of  volume  increase  is  considered  for  a  particular  path. 
For  the  Ironton  Park  path,  snow  entrainment  did  not 
result  in  appreciable  change  in  runout  distance.  The 
additional  momentum  due  to  increased  mass  moving 
down  the  slope  was  apparently  cancelled  by  the 
resistance  of  the  snow  that  was  entrained  during  flow 
over  the  level  runout  zone.  The  primary  benefit  of 
this  study  was  to  set  up  a  mechanism  in  AVALNCH 
for  treating  snow  entrainment. 

There  are  several  notes  of  caution  that  should  be 
mentioned  concerning  the  use  of  program  AVALNCH 

1 .  The  program  predicts  only  the  flow  characteris- 
tics of  the  core  material — not  the  snow  dust 
cloud.  In  cases  where  avalanches  run  beneath 
structures  such  as  bridges,  the  failure  to  ac- 
count for  the  snow  dust  cloud  could  lead  to 
underdesign  of  the  structure. 

2.  Limited  field  observations  indicate  that  long 
runout  distances  are  not  necessarily  associated 
with  extra  thick  slabs  as  indicated  by  the 
program.  In  Colorado,  thick  slabs  are  usually  in 
wind  deposition  areas  and  are  usually  hard 
slabs.  Such  snow  is  strongly  bonded  and  resists 
disintegration  as  it  moves  downslope.  As  a 
result,  flow  height  is  roughly  equal  to  initial  slab 
thickness,  runout  distance  is  not  great,  and  the 
debris  tends  to  mound  up  in  the  runout  zone. 
Fresh  soft  slab,  on  the  other  hand,  is  weakly 
bonded  and  quickly  fluidized  after  fracture. 
This  leads  to  greater  flow  heights,  longer  runout 
distances,  greater  speeds,  and  a  uniform,  shal- 
low layer  of  debris  in  the  runout  zone.  Hence, 
more  field  data  are  needed  to  estabhsh  better 
estimates  of  kinematic  viscosity  and  surface 
friction  for  various  types  of  snow.  The  current 
program  can  handle  this  situation  if  given 
proper  input. 

In  summary,  a  computer  code  has  been  developed 
based  upon  transient  fluid  dynamic  equations  to  model 
snow  avalanche  flow.  In  test  cases  the  program  has 
predicted  results  consistent  with  known  field  observa- 
tions for  a  number  of  avalanche  events.  Capability  to 
modify  the  program  to  treat  special  conditions  such 
as  variable  width  flow,  variable  friction,  snow  entrain- 
ment, and  other  features  has  been  demonstrated. 
Further  refinement  and  application  of  the  program  to 
more  complex  flow  conditions  are  basically  limited 
by  a  lack  of  field  data  on  avalanche  flow  which  can  be 
used  for  verification. 
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Appendix  A 
Input  Data  Format 


The  following  input  data  is  needed  for  program  AVALNCH: 

Card  1:  FORMAT  (20A4) 
Columns  1  through  80:  title  and  identification  information. 

Card  2:  FORMAT  (8F  10.0) 


Columns  1  through  10:  IBAR* 
Columns  11  through  20:  JBAR* 
Columns  21  through  30:  DELX 
Columns  31  through  40:  DELY 
Columns  41  through  50:  NU 
Columns  51  through  60:  FRK** 
Columns  61  through  70:  TWFIN 
Columns  71  through  80:  CWPRT 


—  number  of  cells  in  x  direction. 

—  number  of  cells  in  y  direction. 

—  cell  dimension  in  meters,  x  direction. 

—  cell  dimension  in  meters,  y  direction. 

—  kinematic  viscosity  (m"/s). 

—  surface  friction  coefficient. 

—  avalanche  flow  time  in  seconds. 

—  cycles  between  extended  printouts. 


*Using  version  of  AVALNCH  listed  in  Appendix  B,  IBAR  <  200,  and  JBAR  =  2. 
These  limits  can  be  changed  by  changing  the  dimension  statement  in  the  code. 
**FRK  =  0.0  if  different  values  are  to  be  used  for  different  ceUs  of  the  path. 

**FRK  ^0.0  if  a  constant  value  is  to  be  used  over  the  entire  path  (in  this  case  list 
value  to  be  used). 

Columns  1  through  10 
Columns  11  through  20 
Columns  21  through  30 
Columns  31  through  40 
Columns  41  through  50 
Columns  51  through  60 
Columns  61  through  70 
Columns  71  through  80 


Card  3:  FORMAT  (8F  10.0) 

thickness  of  avalanche  slab  in  cell  1  in  meters, 
thickness  of  avalanche  slab  in  cell  2  in  meters, 
thickness  of  avalanche  slab  in  cell  3  in  meters, 
thickness  of  avalanche  slab  in  cell  4  in  meters, 
thickness  of  avalanche  slab  in  cell  5  in  meters, 
thickness  of  avalanche  slab  in  cell  6  in  meters, 
thickness  of  avalanche  slab  in  cell  7  in  meters, 
thickness  of  avalanche  slab  in  cell  8  in  meters. 


(Continue  on  succeeding  cards  for  IBAR  entrees,  including  zero-thickness  cells.) 

Card  4:  FORMAT  (8F  10.0) 


(Take  these  data  from  longitudinal  profile  of  slope) 


Columns 
Columns 
Columns 
Columns 
Columns 
Columns 
Columns 
Columns 


1  through  10 
11  through  20 
21  through  30 
31  through  40 
41  through  50 
51  through  60 
61  through  70 
71  through  80 
(continue 


change  of  height  from  cell  1  to  cell  2  in 
change  of  height  from  cell  2  to  cell  3  in 
change  of  height  from  ceU  3  to  cell  4  in 
change  of  height  from  ceU  4  to  cell  5  in 
change  of  height  from  cell  5  to  cell  6  in 
change  of  height  from  cell  6  to  cell  7  in 
change  of  height  from  ceU  7  to  cell  8  in 
change  of  height  from  ceU  8  to  cell  9  in 
on  succeeding  cards  for  IBAR  entrees) 


meters, 
meters, 
meters, 
meters, 
meters, 
meters, 
meters, 
meters. 


Card  5:  FORMAT  (8F  10.0) 
(This  set  of  data  required  if  FRK  =  0.0  on  card  4) 


Columns 
Columns 
Columns 
Columns 
Columns 
Columns 
Columns 
Columns 


1  through  10 
11  through  20 
21  through  30 
31  through  40 
41  through  50 
51  through  60 
61  through  70 
71  through  80 


friction  coefficient  in  cell  1. 
friction  coefficient  in  cell  2. 
friction  coefficient  in  cell  3. 
friction  coefficient  in  cell  4. 
friction  coefficient  in  cell'  5. 
friction  coefficient  in  cell  6. 
friction  coefficient  in  ceU  7. 
friction  coefficient  in  cell  8. 


(continue  on  succeeding  cards  for  IBAR  entrees) 


Two  example  data  decks  are  shown  following  the  first  for  friction  specified  in  an 
array,  and  the  second  for  friction  constant  over  the  flow  domain. 
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Appendix  B 

Fortran  Listing  of  AVALNCH— Unit  Width  Version 


1.  C       PROGRAM  AVALNCH. 

2.  ♦ 

3*  C      EQUIVALENT  HORIZONTAL  GRID  OPTION 

4.  ♦ 

5.  C  THIS   PROGRAM  HAS   BEEN   DEVELOPED  TO  NUNERICALLY  SOLVE  (USING 

6.  C  FINITE  DIFFERENCE   TECHNIQUES)   THE  NAVIER  -   STOKES  EQUATIONS 

7.  C  FOR  TRANSIENT  FLUID  FLOW  PROBLEMS.     IT  WILL  BE   USED  TO  SIMULATE 

8.  C  THE  PROBLEM  OF  SNOW  AND  ICE  AVALANCHES  ON  SLOPES  OF  VARYING  SLOPE. 
9«  DIMENSION  U(2  02  t4),V(202,4 ),UN(^02,4),VN(202,4),P(202,4),FR(202  ), 

10.  IX pure   8),H(20  2),HN(202), JTC202  ), NAME(2  0)  ,GX C 20 2 )  ,GY ( 202 )  ,FRCC20 2) 

11.  REAL  NU 

12.  INTEGER  CYCLE 

13.  RE^O    (105,45)  NAME 

14.  WRITE  (108,35) 

15.  WRITE    (108,45)  NAME 

16.  ♦ 

17.  C   ♦  ♦  READ   AND  PRINT    INITIAL   INPUT  DATA 

18.  * 

19.  READ(  105,  25  )   ( XPUT( I ) , 1=1 , 8  ) 

20.  IBflR=XPUT(l>    ;    JBAR=XPUT(2)    i   0ELX=XPUT(3)    ;  0ELY=XPUT(4) 

21.  NU=XPUT(5)      ;      FRK=XPUT(6)      ;      T WF IN=X PUT( 7 )      ;     C WP R T= XPUT ( 8 ) 

22.  WRITE    (  108,  50)   ( XP UT ( I  )  , 1  =  1 , 8  ) 

23.  25  FORMAT(8F10.0) 

24.  35  FO^MAT(lHl) 

25.  45  F0RMAT(?0A4) 

26.  47  F0RHAT(6X'I '7X'J'8X'U'13X'V'liX'P'13X'H'9X,'SUR  CELL') 

27.  48  FORMAT(4X,I3,5X,Ii,4(4X,lPEl0.3),6X,I2) 

28.  49  F0RMAT(2X'CYCLE='I3,2X'ITER='I3,2X'DELT='1PE9.?,2X'TIME='E9.2,2X, 

29.  l'FV0L='E9.2,2X'UMAX='E9.2,2X 'UEDG='E9. ?, 2X'LDEG='I3) 

30.  50  FORMATdH    ,  1 X 'I  B  AR  = 'F  4.  0 ,  2  X  '  JB  AR= 'F3  .  0  ,  2  X  '  DE  LX= 'F6  .2  ,2  X 'OE  L  Y= 'F  5  .  2 

31.  1,2X'NU='F4. 2,2X'FRK='F4.2,2X 'TWFIN='F5.0,2X'CWPRT='F5.0) 

32.  60  FORMAK  8F10.0) 

33.  61   FORMAK  8F10. 3) 

34.  70   FORMAK  IHO,  35X, 'FLOW  HEIGHT') 

35.  71   FORMAKIHO,  25X'ELEVATT  ON  CHANGE  FOR    EACH  CELL') 

36.  72   FORMAKlH0,25X, 'BOUNDARY   FRICTION  COEFFICIENTS') 

37.  73   FORMAK IHO,  25X, 'SLOPE-PARALLEL   GRAVITY  COMPONENT') 

38.  74  FORMAK  IHO,  25  X,  'SLOPE-NORMAL   GRAVITY  COMPONENT') 

39.  75  FORMAK  IHO,  iOX, 'END   OF   INPUT  DATA'//) 

40.  82  FORMAT( 5X, 'PRORLE M  RUNNING  TIME   E X CE E DE D- C A LC U L AT  10 NS  TERMINATED') 

41.  83   FORHATC 5X , ' AV ALANCHE    AT    END   OF   GR ID- C ALC UL AT  I  0 N S  TERMINATED') 

42.  84   FORMAT(SX, 'FLOW   VELOCITY   NEGLIGIBLE-CALCULATION'S  TERMINATED') 

43.  ♦ 

44.  C  *   ♦   COMPUTE   CONSTANT   TERMS   AND    INITIALIZE   NECESSARY  VARIABLES 

45.  ♦ 

46.  lMflX=T'1AR+2    ;    JMAX  =  JBAR*-2 

47.  IM1=IMAX-I      ;  JM1=JMAX-1 

48.  Rrx  =  1.0/DELX 

49.  RDY=1.0/0ELY 

50.  D£LM= DELY/100. 

51.  DELT=1.0 

52.  IM2=IMAX-2 

53.  JM2=JMAX-2 

54.  T=FL&=UE0G1=0. 0 

55.  CYCLE=ITER=INO=LDEG=0 

56.  G=9.8   ;    0MG=1.7    ;    EPSI=.Q01    ;    ALPHA=0.1    ;   GAMMA=0.1    S  DZRO=1.0 

57.  G=9.8   i   0MG=1.7    J    EPSI=.05      ;    ALPHA=0.1    ;   GAMMfl=0.1    i  DZRO=1.0 

58.  BETA=   0MG/(  2.*DELT♦(RDX♦♦2-^RDY*»2)) 

59.  ICPRT=INT(CWPRT) 

60.  IF(ICPRT.EQ.l)  ICPRT=2 

61.  DO  100  I=1,IMAX 
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62. 

H(I)  =  HNCI)= JT( I)  =  GXCI)=GYC I)=FRC(  I)  =  0.0 

63. 

00   100  J=ltJMAX 

64. 

100 

UCI» J)=V(I, J)=UN(I, J)=VN(I,J)=PC ItJ)=0.0 

65. 

t 

66. 

C 

♦  * 

SPECIAL   INPUT  DATA 

67. 

♦ 

68. 

READC 10  5,60)(HCI)tI=2,IMl) 

69. 

REAOC 10  5,60)(HNCI)t 1  =  2  ,T«1 ) 

70. 

IFCFRK.GT.O.O)   GO   TO  120 

71. 

REA0(105,60)  (FRCCI)tI=2,IMl) 

72. 

GO  TO  130 

73. 

120 

00  125  1=2, IMl 

74. 

125 

FRCCI  )=FRK 

75. 

130 

CONTINUE 

76. 

00  15  0  I=2tl«l 

77. 

FR(I)=FRC<I) 

78. 

SP=HN(I)/0ELX 

79. 

CPsSQRTC1.0-SP*SP) 

80. 

GX(I)=G*SP 

81. 

150 

GY(I)=-G«CP 

82. 

WRITE<106.70; 

83. 

WRITE  CI 08f61)(H(I)»I»2, IMl) 

84. 

URITEC108,71) 

85. 

WRITE(108f 61)(HNCI)tI^2«INl) 

86. 

WRITE(108,72) 

87. 

tfRITEC108,61)  (FRCCI), 1=2, IMl) 

88. 

WRITEC108,73) 

89. 

URITEC108,61)CGXCI),I=2,IM1) 

90. 

WRITEC108,74) 

91. 

WRITEC108,61)(GYCI)f I=2,1M1) 

92. 

WRITEC108,75) 

93. 

DO  240  1=2, IHl 

94. 

JTCI)=INTCHCI )*R0Y+l.E-6)*2 

95. 

IFCJTCD.GT.JMl)  JTCI)=JM1 

96. 

240 

HNC  I)=0  .0 

97. 

HC1)=  HC2) 

98. 

H(IMAX)=  HCIMl) 

99. 

JTC1)=  JT(2) 

100. 

JTCI«AX)=  JTCIf^l) 

101. 

102. 

c 

«  « 

CALCULATE   HYDROSTATIC  PRESSURE 

103. 

104. 

DO  280   1=2, IMl 

105. 

JT1=JTC I) 

106. 

DO  280  J=2,JT1 

107. 

280 

PCI, J)=-GY( I)*CHCI)-(FL0AT(J)-1.5)*DELY) 

108. 

ASSIGN  4280   TO  XRET 

109. 

GO  TO  2000 

110. 

* 

111. 

c 

*  « 

START  CYCLE 

112. 

* 

113. 

1000 

CONTINUE 

114. 

ITER=0 

115. 

FLG=1. 

116. 

ASSIGN   3000   TO  KRET 

117. 

118. 

c 

COMPUTE   TEMPORARY  U   AND  V 

119. 

120. 

00   1100   1=2  , IMl 

121. 

JT1=JTC  I) 

12?. 

on   1100  J=2,JT1 

1?3. 

FUX=(CUN(I,J)+UN(I+l,J))*(UN(I,J)+UN(I*l,J))*ALPHA*ABS(UN(ItJ)*UN( 

124. 

ll^l, J))*(UN(I, J)-UN(I*1, J))-(UN(I-1, J)*UNCI,J))*CUNfl-l,J)*U\(r ,J) 

125. 

2)-ALPHA*ABSCUNCI-l,J)*U.N(I,j;)«(UNCl-ltJ)-UN(I,J>J)/(4.*0ELXi 

126. 

FUr=( CVNC I, J)*VNCI*1, J))*CUN(I , J)*UN(I , J*l)) 
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127, 

l*ALPHA*ABS(VN(I,J)+VN(I+l,J))»(UN(ItJ)-UN(ItJ*l)) 

128. 

2-(VN(  I,J-l)*VN(I+l,J-l))*CUN(I,J-l)+UN(IfJ)) 

129. 

3-ALPHA*ABS(VN(I,J-l)+VS(I  +  lfJ-l))*CUN(I,J-l)-UN(I,J)))/(4.*0ELY') 

130. 

FVX  =  C(UN(I,J>*UN(r,J*l>;*CVNCr»J)+VN(.I*l,J))*ALPHA*ABS(UNCIfJ)  +  UN( 

131  . 

lI»J+l))*(VN(IfJ)-VN(I+ltJ))-(UN(I-l,J)+UNCI-l»J+l))*CVNCI-l»J)*VNC 

132. 

*?I,J))-ALPHA*ARSCUNCI-1»J)*UN(I-1,J+1))*(VN(I-1,J)-VN(I,J)))/(4.*0E 

133. 

3LX:) 

134. 

FVY=((VN(I,J><^VNClfJ  +  l))«CVN(I,  J:)+VNC  I  ,  J  +  1  )  )  ♦A  L  PH  A*  AB  S(  VN(  I,  J)  +  VN 

135. 

l(ItJ*l))*CVN(I,J)-VN( I , J+l))-(VNCIf J-1)*VN( I, J))*CVN(If J- U^VNCIf J 

136. 

2))-ALPHA*ABSCVNCIfJ-l)  +  VN(I f  J ) )♦ ( V N ( I » J- 1 ) - VN C I tJ)))/C4.*0ELY) 

137. 

VISX=   NU*((UNCI*lf J)-2.*UN(I,J)*UN<I-1  ,J))/0ELX**2* 

138. 

1                      CUNCI»J*l)-2.*UN(I,J)*UN(It J-1))/0ELY**2) 

139. 

VISY=  NU*((VN(I*ltJ)-2.*VN(I,J)*VN(I-ltJ))/DELX**2» 

140. 

1                      ( VNCI, J+li-2.»VN(I , J)*VN(It J-1))/0ELY**2) 

141. 

U(I,J)=   UN( If J)*OELT*((P(Tf J)-P(I*lf J))#RDX   ♦  GX(1 )-FUX-FUY*VIS X) 

142. 

1100 

VCI»J)=  VNCI,  J)*OELT*(<.P(I,J)-PCIt  J*1>3*RDY   ♦   GY(  I  )-FVX-FVY  +  V  ISY) 

143. 

144. 

c  ♦  * 

SET   BOUNDARY  CONDITIONS 

145. 

♦ 

146. 

2000 

CONTINUE 

147. 

HN(1)=  mi2) 

148. 

HN( IM AX)=HN(IM1) 

149. 

JTC1)=JTC2) 

150. 

JT(I«AX)=JT(IM1 ) 

151  . 

C       LEFT  WALL  RIGID   AND   SLIP  FREE. 

152. 

C       RIGHT   WALL  CONTINUOUS  OUTFLOW. 

153. 

DO  2200  J^lfJMAX 

154. 

UClt J>=0.0 

155. 

V(l ,J)=VC2, J3 

156. 

IFCITER.GT.O)  GO  TO  2200 

157. 

U(IN1  ,J)=U(T«2, J) 

158. 

VCinAXt J)=VCINlt J) 

159. 

2200 

CONTINUE 

160. 

C       TOP  WALL  CONTINUOUS  OUTFLOW. 

161  • 

C       BOTTOM  WALL  RIGID  -  WITH  FRICTION 

162. 

00  2500  I=l9lMAX 

163. 

IFCITER.GT.O)  GO  TO  2400 

164. 

VCIff JHl)»VCIf JM2) 

165. 

U(ItJHAX)«UCIf JNl) 

166. 

2400 

V(lfl)«0.0 

167. 

2500 

U(Ifl>«U(I»2)*C1.0-2.0*FRCCI>) 

168. 

169. 

C  ♦  ♦ 

FREE    SURFACE  BOUNDARY  CONDITIONS 

170* 

171. 

DO  2650  I:=:29lMl 

172. 

JT1=JT(I) 

173. 

IF(  JTd^D.LT.JTd))  U(IfJTl)=U(  I»JT1-1) 

174. 

V(I,JT1)=  VCIt  JTl-l)-DELY*RDX*(UCItJTn-UCI-l,  JTD) 

175. 

2650 

UCI.JT1*1)=  UCI.JTI) 

176. 

GO  TO  KRET,C3000t4280> 

177. 

3000 

CONTINUE 

178. 

♦ 

179. 

c  ♦  ♦ 

HAS  CONVERGENCE   BEEN  REACHED 

180. 

* 

181. 

IFCFLG.EO.O .)   GO  TO  4000 

162. 

ITER=ITER*1 

183. 

IF(ITER.LT.500)  GO  TO  3050 

184. 

IF(CYCLE.LT.IO)  GO  TO  4000 

185. 

T  =  1.E*10 

186. 

GO  TO  4000 

187. 

3050 

FLG=0.0 

188. 

189. 

c  ♦  ♦ 

COMPUTE  UPDATED  CELL  PRESSURE   AND  VELOCITIES 

190. 

* 

191. 

JB1  =  2 
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192.  00  J500   1=2, I«l 

193.  JTi=JT(T) 

194.  DO  3500  J=2,JT1 

195.  IFC JTl.EO. JBl)   GO   TO  3G60 

196.  IFCJ.NE.JBl   .AND.    J.NE.JTl)   GO  TO  3200 

197.  IF(J.EQ.JTl)  GO   TO  3100 

198.  GO  TO  32C0 

199.  3060  CONTINUE 

200.  F=V(I  ,J)>DELY*RDX*(UCI,J)-UC I-l, J)) 

2  01.  0F0P=DELT*R0Y*(1. 0*2. 0*DELY*»2*R0X**2) 

202.  DELP1=-F/DFDP 

203.  3100   ETft=DELY/CHNCIJ-CFLOATC JTl>-2. 5>*DELY) 
20  4.  D£LP=(1.0-ETA)*P(I,JT1-1)-P(I,JT1) 

205.  IF( JB 1. EQ.JTi  )   DELP  =  a.5*(nELP*0eLPl) 

206.  GO   TO  3300 

207.  3200   D=RDX*CU(T,  J>-UCI-l,J)>  +  RDY»(V(r,  J^-VCIjJ-D) 
20P.  IF( 4B SCD/DZRO).GE.EPSI)  FLG=i.O 

209.  n£LP=  -PETA*0 

210.  3300   P( I , J )=P( I, J)4DELP 

211.  UCT,J^=UCI, J>+DELT*R0X*0ELP 

212.  U( I-l  ,J)  =  U(  I-l , J)-DELT*ROX*DELP 

213.  V(I,J)=V( I, J)*DELT*RDY*DELP 

214.  3500   VCr , J-1)=V(T, J-1)-DELT*R0Y*DELP 

215.  GO   TO  2000 

216.  4000  CONTINUE 

217.  ♦ 

218.  C   *  ♦  COHPUTE   NEW   POSITION   FOR   TOP  SURFACE 

219.  ♦ 

220.  DO  4100   1  =  2, IMi 

221.  JT1  =  JTCI) 

?22.  HV=  RDY*(HN(I)-FL0AT(JT1-2)*0ELY) 

223.  UAV=   0.5*(UCI-1, JTl)   ♦  UCI»JT1)^ 

224.  HCI)=HNCI  )*FVOLl/FVOL+DELT*(.HV*V(I,JTi;  +  Cl  .0-HV;*VCI,  JTl 

225.  l-l)-0  .5*R0X*CUAV*HNCI-H)*GAMMA*aos(UAV)*(HN(  D-HNd  +  D) 

226.  2       -U  AV*HN(  I-l)-GA^'^^A♦  AeS(UAV)«(HNCI-l)-HNCI)))) 

227.  4100  CONTINUE 

228.  ♦ 

229.  C  ♦   ♦   CALCULATE   CELL   IN  WHICH  SURFACE    IS  LOCATED   AND  UPDATE  ARRAY 

230.  ♦ 

231.  00  4250   1=2, IMI 

232.  IF(H(I).LT.DELM)  H(I)=0.0 

2  33.  JTCI>=INTCHCI)*R0Y*1.0E-6>*2 

234.  IF(JTCI).GT.JHl)  JT(I)=JH1 

235.  4250  CONTINUE 

236.  ASSIGN  4280  TO  KRET 

237.  GO  TO  2000 

238.  4280  CONTINUE 

239.  ♦ 

240.  C  ♦  ♦  CALCULATE  TOTAL  FLUID  VOLUME 

241.  ♦ 

242.  FVOL*0. 0 

243.  DO  4300  I=2flMl 

244.  4300  FVOL=FVOL*H(I)*OELX 

245.  IF(CTCLE.EQ.O)  FV0L1=FV0L 

246.  ♦ 

247.  C  ♦  •  FIND  LEADING   AND  TRAILING  EDGES  OF   AVALANCHE  AND  LD.   EDGE  VELOCITY 

248.  • 

249.  LOEGl^DEG 

250.  I«IH2 

251.  4400  IF(H(  D.GT.DELM)  GO  TO  4500 

252.  1=1-1 

253.  GO  TO  4400 

254.  4500  LDEG«I 

255.  1*2 

256.  4600  IF(HCI).GT.DELM)  GO  TO  4700 
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257.  1=1*1 

258.  GO  TO  4600 

259.  4700  KTEG*I 

260.  IF(LDEG.EQ.LOEGl)  GO  TO  4800 

261.  IFCCYCLE.GT.O)  UEOG=DELX/TC 

262.  IF( CYCLE. EG. 0)  UEDG=5.0 

263.  TC=OELT 

264.  INFL0=1 

265.  IF(UEDG.GT.UEOGl)  UEDG1=UEDG 

266.  GO  TO  4910 

267.  4800  TC=TC*DELT 

268.  INFL0=INFL0*1 

269.  ♦ 

270.  C  ♦  ♦   ADVANCE  UtV.H  ARRAYS. 

271.  * 

272.  4910  UMAX=VHAX=0.0 

273.  00  4900  I=1»IMAX 

274.  00  4900  J=1,JNAX 

275.  IFCABSCUCI, J^J.GT.1.0E*04)  UCI,JJ=0.0 

276.  UNC It J)=U( I 

277.  IF( ABSCVCIt J)).GT.1.0E*04)  VCItJ)=0.0 

278.  VN(I, J)=VCI t J) 

279.  IFC AR SCPCIt J)).LT.l .0E-16>  P(I,J)=0.0 

280.  4900  HNCI)=H(I) 

281.  DO  4950  I=KTEG,L0EG 

282.  DO  495U  J=2»JM1 

283.  UT  =  ABS(UN(I  ,J)) 

284.  VT=ABSCVN(I , JJ^ 

285.  IFCUT .GT.UMAX)  UMAX=UT 

286.  4950   TFC VT .GT. VMAX)  VMAX=VT 

287.  ♦ 

288.  C   ♦   ♦  LIST  VELnCITY,    PRESSURE,   AND   SURFACE  POSITION 

289.  ♦ 

290.  5000   WRITE(1  08,49)   C  YCLE  »  I  T  E  R  ,  OE  L  T  ,  T  ,  FVOL .  U   AX  ,  UE  DG  ,  LDEG 

291.  IFCCYCLE.EQ. ICPRT)  GO   TO  5030 

292.  IFCCYCLE. NE .ICPRT)   GO   TO  600C 

293.  5030  ICPRT=   ICPRT   ♦  INT(CWPRT) 

294.  5060  CONTINUE 

295.  WRITE (108,47) 

296.  DO  5250  I=1,IWAX 

297.  JT1=  JTCI) 

298.  JT2=JT1*1 

299.  DO  5250  J=1,JT2 

300.  WRlTE(i08,48)    I , J ,UCI , J) , V( I , J ) , P C I , J  J , HCI ) , J Tl 

301.  5250  CONTINUE 

302.  IFCIND.EQ.2)  GO  TO  6520 

303.  IF(IND.EQ.3)  GO  TO  6530 

304.  IF(IN0.EQ.4)  GO   TO  6540 

305.  ♦ 

306.  C  ♦  ♦  RECOMPUTE  CONTROL  PARAMETERS. 

307.  » 

308.  6000  IFCCYCLE. EO.O)  GO  TO  6300 

309.  DTX=OELX/UNAX 

310.  DTY=OELY/VMAX 

311.  OELT= A«IN1(OTX,DTY)/3.0 

312.  IFCITER.lt. 10)  DELT=1 . 5*DELT 

313.  IFCNU-1.E-6.lt. 0.0)  GO  TO  6300 

314.  0ET=CDELX*DELY)»*2/C2.*NU*CDELX«*2*DELY**2>) 

315.  IFCDELT.lt. DET)  GO  TO  6300 

316.  DELT=0.9*OET 

317.  6300  T=T*0ELT 

318.  IFCCYCLE. EQ.O)  GO  TO  6400 

319.  0AX»UHAX«DELT/0ELX 

320.  DAY=V«AX*OELT/DELY 

321.  ALPHA=1.35«AMAX1C0AX,0AY) 
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322. 

IFCALPHA.GT.l.G)  ALPHA=0.95 

323. 

GAHMA^ALPHA 

324. 

BETA*0HG/(2.*0ELT*(R0X**2*R0Y**2)) 

325. 

* 

326. 

C  ♦  ♦ 

TEST  FOR  PROGRAM  TERMINATION. 

327. 

328. 

6400 

1F( T.GT.TWF  IN)  lNn=^ 

329. 

IFCHCIN2J  .GT.OFLM)    IND  =  3 

330. 

IF(UED&.LT.0.05*UEDG1)  IN0=4 

331. 

IFCINFLn.EQ.50)  IND=4 

332. 

IFCIND.GT.i)  GO  TQ  6500 

333. 

IFCCYCLE.lt  .3)  GO   TO  644C 

334. 

AA  =  1.0*20.0*tXP(-l.^:5*UEDG; 

335. 

00  64  30   I  =  2.I'^1 

336. 

6430 

FRC(I)=Afl*FRCT) 

337. 

6440 

CYCLE=CYCLEtl 

338. 

GO  TO  lOCO 

339. 

6500 

T=T-DELT 

340. 

GO  TO   5  06  0 

341. 

6520 

WRITE (108, 82) 

342. 

GO   TO  6600 

343. 

6530 

WRITEC108,83) 

344. 

GO  TO  6600 

345. 

6540 

WRITEC108,84) 

346. 

6600 

STOP 

347. 

ENO 
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Appendix  C 
Listing  of  Output  Data  from  AVALNCH 
for  Ironton  Park  Avalanche 


IROMOK    PSPK    aVlLANCHE    PATH   -   HORlZOMfiL    GRID  OPTIO\. 
ie*R  =  100.      jeAR-   2.      OELX=   IC.OO     DELY=   2.00      NU=  .5C      FRK=    .50      TWFU=  750  .     CWPRT  =  lOCC. 


FLCW  HEIGHT 


2.000 

2.000 

2.  OOC 

2.  CCC 

i^.OOO 

1.500 

l.COC 

.500 

.OCO 

.000 

.  OUC 

.OOC 

.000 

.000 

.000 

.000 

.OCO 

.000 

.  000 

•  OCC 

.000 

.coo 

.000 

.COC 

.000 

.000 

.COD 

.000 

.000 

.000 

.000 

.000 

.  occ 

.  000 

.  OOC 

.CCC 

.000 

.coo 

.000 

.OOC 

.OCO 

.000 

.  000 

.OCO 

.000 

.COO 

.000 

.000 

.OCO 

.  000 

.  OOC 

.000 

.000 

.000 

.000 

.000 

.000 

.000 

.OOC 

.OCO 

.000 

.000 

.000 

.000 

.ooc 

.  000 

.  OOC 

•  OCO 

.000 

.COO 

.  OOC 

.000 

.000 

.000 

.OOC 

.OCO 

.000 

.000 

.000 

.000 

n  fi  n 

•  u  \J  U 

n  n  n 

.000 

.000 

u  n  n 

.  u  u  u 

n  n  f) 

fi  n  n 

r  r  f 

.  U  U  1/ 

.000 

.000 

.000 

.OCC 

.000 

.000 

.000 

.000 

.000 

.000 

.  000 

•  OCO 

ELfcVftTION  CHANGE 

FOR  EACb 

CELL 

6. 500 

6.200 

6.200 

6.200 

6.500 

6.2C0 

6.200 

6.200 

6.500 

6.500 

6.200 

5.800 

b.OOO 

5.800 

5.600 

5.600 

5.200 

4.600 

4.  ?0C 

4.200 

4.600 

5.200 

4.600 

4.800 

5.000 

5.200 

5.400 

5.600 

5.800 

6.200 

6.200 

6.200 

6.200 

6.500 

6.  800 

7. 400 

7.200 

7.200 

7.100 

7.000 

7.  OCC 

7.000 

6.600 

6.  2C0 

6.200 

6.000 

6.000 

5.600 

5.2CC 

5.^00 

5.000 

4.600 

4. 200 

4.200 

4.000 

4.000 

3.6C0 

3  .  6C0 

3.400 

3. 400 

3.000 

2.800 

2.500 

2.200 

1.8CC 

1.600 

1.400 

1.200 

.900 

.700 

.500 

.4C0 

.300 

.100 

.000 

.OCC 

.000 

.COO 

.000 

.000 

.OCC 

.000 

.  000 

.OOC 

.000 

.CCC 

.000 

.000 

.000 

.000 

.000 

.000 

.000 

.COO 

.OOC 

.000 

.000 

.000 

.000 

.OCO 

BCUNOftPY  FRICTION 

COEFFICIENTS 

.500 

.500 

.500 

.5C0 

.500 

.bOO 

.500 

.500 

.SCO 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.5C0 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.  500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.5,0  0 

.  500 

.500 

.500 

.500 

.500 

.500 

.500 

.5  00 

.500 

.5C0 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

•  ^  t  u 

•  7  u  u 

.500 

.5C0 

.500 

.500 

•  7  u  U 

i\  n  r 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

.500 

SLOPE 

-PARALLEL  GRAVITY  COHPOKENT 

6.370 

6.076 

6.  076 

6.076 

6.370 

6.076 

6.076 

6.076 

6.370 

6.  370 

6.  076 

5.684 

5.880 

5.684 

5.488 

5.488 

5.056 

4.508 

4.116 

4.116 

4.508 

5.096 

4.508 

4,704 

4.900 

5.096 

5.  292 

5.488 

5.684 

6.076 

6.076 

6.076 

6.076 

6.370 

6.664 

7.056 

7  .056 

6.958 

6.860 

6.860 

6.860 

6.468 

6.076 

6.076 

5.880 

5.880 

5.488 

5.096 

5.096 

4.900 
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3.920 

3.920 

3. 528 
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2.940 
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.882 
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.  000 

.000 

.000 
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.000 
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.  000 

.000 
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.000 

.COC 

.000 

.000 

.000 

.000 

.000 

.000 

.000 

.OOC 
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.000 
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.  OOC 

.OCO 

SLOPE 

-NORMAL   GRAVITY  CO^•PO^E^T 

-7.447 

-7.689 

-7.689 

-7.689 

-7.447 

-7.689 

-7.689 

-7.689 

-7.447 

-7.447 

-7.689 

-7.983 

-7.840 

-7.983 

-8.119 

-8.119 

-8.371 

-8.702 

-8.  894 

-8.894 

-8.702 

-8.371 
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-7.689 
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-9.560 
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-9.760 
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-9.792 
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-9.800 

-9.  800 

-9.800 

-9.800 

-9.800 

-9.800 

-9.800 

-9.8CC 

-9.800 

-9.80C 

-9.8C0 

END   OF    INPUT  RATA 
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01 
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01 
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01 
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nE  LT  = 
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01 

yifE  = 
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01 

FVCL  = 

1 .  30E 

02 
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UECG* 
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01 
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99 

ITER  = 
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01 
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1 

-4.635F- 

01 

.OOOE 

00 

.OOOE 

00 

2.204E 

00 

3 

94 

2 

3.271E 

00 

-4.633E- 

02 

1.194E 

01 

2.204E 

00 

3 

94 

3 

5.0  16E 

00 

-1.620E- 

01 

7.596E 

00 

2.204E 

00 

3 

94 

4 

5.0  16E 

00 

.OOOE 

00 

.OOOE 

00 

2.204E 

00 

3 
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95 

1 

-4.838E-01 

.OOOE  00 

.OOOE 

00 

2.063E 

00 

3 

95 

2 

3.415E  00 

1  .834E-02 

1.026E 

01 

2.063E 

00 

J 

95 

3 

5.027E  00 

1.620E-02 

-9.468E 

00 

2.063E 

00 

3 

95 

4 

!>.027E  00 

.OOOE  00 

.OOOE 

00 

2.063E 

00 

3 

96 

1 

-7.060E-01 

.COGE  00 

.OOOE 

00 

2.107E 

00 

3 

96 

2 

4.983E  00 

-3.261E-01 

1  .163E 

01 

2.107E 

00 

3 

96 

3 

4.983E  00 

-3.173E-C1 

-7.458E 

00 

2.107E 

00 

3 

96 

4 

4.983E  00 

.OOOE  00 

.OOOE 

00 

2.107E 

00 

3 

97 

1 

-2.1 86E-01 

.OOOE  00 

.OOOE 

00 

1 .561E 

00 

2 

97 

2 

1.543f  00 

6.8e0F-01 

-4. 201F- 

03 

1  .561E 

00 

2 

97 

3 

1.543F  00 

-J.184E  CO 

.OOOE 

CO 

1.561E 

00 

2 

98 

1 

-7.aS3F-03 

.OOOE  00 

.OOOE 

00 

2. 8915- 

02 

2 

98 

2 

';.5  75£-02 

2.974E-C1 

5.575E- 

05 

2.891E- 

C2 

2 

98 

3 

5  .5  75F-02 

-4.429E  00 

.OOOE 

00 

7.891E- 

02 

2 

99 

1 

-I  .7576-05 

.CCCE  CC 

.OCOE 

00 

.OOOE 

CC 

2 

99 

2 

1  .240F-04 

1. 112F-02 

8.604F- 

07 

.OOOE 

00 

2 

99 

3 

1  .240E-04 

-5.141E  CC 

.  COCF 

00 

.OOOE 

00 

2 

100 

1 

S  .890  E-08 

.OOOE  CC 

.  COOE 

CC 

.OOOE 

00 

2 

1  00 

^ 

-4.1 57E-07 

2 .488E-C5 

-8.514E- 

cs 

.OOOE 

00 

2 

100 

3 

-4. 1 57E-07 

-5.156E  00 

•  OOOE 

oc 

.OOOE 

00 

2 

101 

1 

b . 8  eH£- C8 

.OOOE  00 

.OOOE 

oc 

.OCCE 

00 

2 

101 

2 

-4.1 56C-07 

-1 .811 F-11 

-1  .031 E- 

U9 

.OCOE 

CO 

2 

101 

3 

-4.156E-07 

-5.182F  CO 

.OOOE 

00 

.OOOE 

CC 

2 

IQZ 

1 

.OOUF  00 

.OOOE  00 

.OOOE 

00 

.OOOE 

00 

2 

102 

2 

.OOOF  00 

-5.182E  OG 

.COOE 

00 

.OCOE 

00 

2 

102 

J 

.UUJF  CO 

-5.182E  00 

.COCE 

00 

.OOOE 

00 

2 
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Appendix  D 
Listing  of  Program  AVALNCH— 
Variable  Width  Option 


1.  C       PROGRAM  AVALNCH. 

2.  ♦ 

3.  C       EQUIVALENT  HORIZONTAL  GRID  OPTION   :    VARIABLE   WIDTH  OPTION 

4.  * 

5.  C  THIS   PROGRAM  HAS  BEEN  DEVELOPED   TO  NU^^ERICALLY   SOLVE  CUSING 

6.  C  FINITE   DIFFERENCE   TECHNIQUES)  THE  NAVIER   -   STOKES  EQUATIONS 

7.  C  FOR   TRANSIENT   FLUID  FLOW  PROBLEMS.      IT   WILL   BE   USED  TO  SIMULATE 

8.  C  THE  PROBLEM  OF  SNOW  AND  ICE   AVALANCHES  ON   SLOPES  OF  VARYING  SLOPE. 

9.  DIMENSION  UC202,7),V(202,7),UN(20  2,7),VN(20  2,7),PC202,7),FR<2  02), 

10.  IXPJTC    8),HC202),HN(2  02)f JT(202)fNAME(2C),GXC2  0Z) ♦GYC202),FRCC23  2)  , 

11.  2W(202) 

12.  REAL  NU 

13.  INTEGER  CYCLE 

lA.  READ   C1C3,45)  NAME 

15.  WRITE  (108.35) 

16.  WRITE   C10Bt45)  NA^E 

17.  * 

18.  C  *   *   READ    AND   PRINT    INITIAL    INPUT  DATA 

19.  * 

20.  READC105,?5)   ( XPUT C I) , I =1 f 8 ) 

21.  IBAR=XPUT(1)   ;    JBAR=XPUT(2)    ;    0ELX=XPUTC3)    J  0ELY=XPUTC4) 

22.  NU=XPUT(5)      J      FRK=XPUTC6)      i     TWFIN=X PUTi 7 )      ;     CWPRT=XPUT< 8) 

23.  WRITE   C108,50)   CXPUTC I) , 1=1 , 8) 

24.  25  FORMAT(8F10.0) 

25.  35  FORMAT(lHl) 

26.  45  FD^MftT(20A4) 

27.  47  F0RMATC6X'I'7X'J'8X'U'13X'V'13X'P'13X'H'9X,'SUR  CELL') 

28.  48  FORMATC 4X, I3t5X, I3f 4(4X ,1PE1 0. 3) »6X f 12) 

29.  49   F0-iMAT(2X'CYCLE  =  'I3,ZX'IT£R"='I3,2l<'DELT=  •'1PE9.2,2X'TIME  =  'E9.2,2X, 

30.  l'FV0L='E9.2,2X'UMAX='E9.2,2X'UE0G='£9.^,2X'L0EG='I3) 

31.  50  FORMATCIH   , IX 'I B AR= 'F4 . 0 , 2 X ' J3 AR ='F 3. 0 , 2X ' DELX= 'F6 . 2 , 2 X 'DEL Y= 'F 5. 2 

32.  1,2X 'NU= 'F4.2|2X'FRK='F4.2,2X'TWFIN='F5. 0, 2X 'CWPRT-  'F5  .0) 

33.  60   FORMATC 8Fi0.0) 

34.  61  FORMATC8F10.3) 

35.  70   FORMATClHOt 35X, 'FLOW  HEIGHT') 

36.  71    FORMAT(lH0,25X'ELEVATinN  CHANGE  FOR    EACH  CELL') 

37.  72   FORMATC IHO, 25X, 'BOUNDARY  FRICTION  COEFFICIENTS') 

38.  73  FORMATC IHO. 25X, 'SLOPE-PARALLEL   GRAVITY  COMPONENT') 

39.  74   FORMATC 1H0,25X, 'SLOPE-NORMAL  GRAVITY  COMPONENT') 

40.  75   FORMATCIHO, 30 X, 'END   OF   INPUT  DATA'//) 

41.  76   FORMA TC IHO ,30X  , 'PER   UNIT   FLOW  WIDTH') 

42.  82   FORMATC 5X , 'PROBLEM   RUNNING  TIME   E XCE E DE D-C A LCU L AT 1 0 NS  TERMINATED') 

43.  83   F0RMATC5X , 'AVALANCHE    AT   END   OF   GRID-CALCULATIONS  TERMINATED') 

44.  84   F0RMATC5X,'FL0W  VELOCITY  NEGLIGIBLE-CALCULATIONS  TERMINATED') 

45.  * 

46.  C   *  *   COMPUTE    CONSTANT   TERMS   AND   INITIALIZE   NECESSARY  VARIABLES 

47.  * 

48.  IMAX=IBAR*2    ;  JMAX=J8AR>2 

49.  IM1=IMAX-1      ;  JMl=JMAX-i 

50.  RDX=1.0/DELX 

51.  RDY=1.0/DELY 

52.  OELM=DELY/100. 

53.  DELT=1.0 

54.  lM2=lMAX-2 

55.  JM2=JMAX-2 

56.  T=FLG=UEDG1=0.0 

57.  CYCLE=ITER=IND=LOEG=0 

58.  &=9.8   5   0MG=1.7   ;    EPSI=.001   i    ALPHA=0.1   J    GAMMA=0.1   ;  DZRO=1.0 
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59.  3ETA=   0»^G/(  2,*DELT*(R0X**2+RDY*»2)) 

60.  ICPRT=INTCCWPRT) 

61.  IFCICPRT.EO.l)  ICPRT=2 

62.  DO  100  I=1,IHAX 

63.  HCI)=HN(I)=JTC I)=GX(I)=GyCI)=FRC(I)=WCI)=0.0 

64.  DO   100  J=ltJMAX 

65.  100   UCIf J)=VCIf J)=UNCI,J)=VNCI,J)=PCI,J)=0.0 

66.  ♦ 

67.  C   *  ♦   SPECIAL   INPUT  DATA 

68.  * 

69.  REAOC105,60>CHCI),I=2,IW1) 

70.  READC105t603CHN(I),I=2tIHl) 

71.  REA0C105t60)CWCI)tI=2,lHl) 

72.  IFCFRK.GT.0.0)  GO  TO  120 

73.  REA0C105,60)  CFRCC I )f I =2, I  Ml ) 

74.  GO  TO  130 

75.  120  DO  125  1=2, IMl 

76.  125  FRC(I)=FRK 

77.  130  CONTINUE 

78.  DO  150  1=2, INI 

79.  frci)=frcci) 

80.  sp=hnci)/delx 

81.  cp=sqrtci.o-sp*sp:) 

82.  GXCI>=G*SP 

83.  150  GYCI)=-G*CP 

84.  WR[TEC108,70) 

85.  WRITEC108,61>CHCI),I=2,IN1) 

86.  WRITEC108,71) 

87.  WRITE(108,61)(HN(I), I=2,IM1) 

88.  MRITE(108,72) 

89.  MRITEC108,61)  (FRCCI), 1=2, INI) 

90.  WRITEC108,73) 

91.  WRITEC108,61)(&X(I),I=2,I>^1) 

92.  WRITEC108,74) 

93.  WRITE(108,61)(GYCI),I=2,IM1) 

94.  WRITEC108,76) 

95.  WRirEClU8,61)CWCI),I=2,INl) 

96.  WRITE(108,75) 

97.  DO  240  1=2, INI 

98.  JTCI>  =  rNTCH(U*R0Y*l.E-6J  +  2 

99.  IFC  JTCD.GT.JNl)  JTCI)=JN1 

100.  240  HN(I)=0.0 

101.  H(i>=  HC2) 

102.  W(1)=W(2) 

103.  HtINAXJ=  H(INl) 

104.  JT(1)=  JT(2) 

105.  JTCIMAX)=  JTCINl) 

106.  ♦ 

107.  C  ♦   ♦  CALCULATE  HYDROSTATIC  PRESSURE 

108.  » 

109.  DO  280   1=2, IMI 

110.  JT1=JTCI) 

111.  00   280  J=2,JTi 

112.  280  PCI,J)=-GYCI)*CHC n-CFL0ATCJ;-1.5)*0ELY> 

113.  ASSIGN  4100   TO  KRET 

114.  GO  TO  2000 

115.  ♦ 

116.  C   ♦  ♦   START  CYCLE 

117.  » 

118.  1000  CONTINUE 

119.  ITER=0 

120.  FL&=1. 

121.  ASSIGN   3000    TO  KRFT 

122.  ♦ 
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171 
X  £.J  • 

r  *  * 

rriMPlITF    TFMPnR&RY    U    AND  V 

* 

IPS- 

DO    110 D  T=2.TM1 

L/C     XXWV/     x  ^9i.rix 

A  CO  • 

nn  1  1  f)fi    {=.?  -  IT  1 

ICO* 

A  t  7  • 

Xi'X^«Jyv*vUMVit»/y    Liiivx^Xf«jyy    vvjnivx    Xf«jy^»ji>iVXf«jyy*vunivx  xjjy^^vjpq^^jjy 

1  ^  n 

cv    Mi-.rn'A*MDOVvjiivx    xf*i/"wnivxfwy-/^vui'ivx     xf^y  W''ivxfvyyy#v^.*ucL.Ay 

A  J  X  • 

FlJY  =  CfVNCI«J)*VNfI*l«J^)*fUNfT«J^*UNCI«J*l"i> 

1    \JI"~V  V*  *"V  Xy            •  »  ''V  X  "X  fwyy^V^'*V^  Fwx*  v/'^VX  Jw*   X  y  / 

1  ^? 
X  JC  • 

l«-flLPHa*flBSfVNfI.J)*VNd-H.J')>»(UNCT.J)-UN(I«J*ll) 

X"Wwri«H*  Muov  wmv  X  Pw           ''X  X  "x  J  «/  J  y^Vw<"Vx  J  wy    \Ji^vx^»/  '  X  yy 

A  J  J  • 

?-(VN(  I«J-1)+VN(  I*1«J-1^J*(UNCI«J-1)*UN(I.J>) 

£-     v»'»vxjw     xy'«i«vx*xyw     xyy^xW'"vxyw  xy*v#*^vxywyy 

X  J  T  • 

^-flLPHA*ARSfVN(I«J-l>  +  VN(I+l«J-l))*CUN(I«J-l')-UN(T.J^^^/f4-*DFLY") 

1 

FVX=ffllWrT-  l^*lJNfT-J+l^')*fVNfr-J'>  +  VNfT  +  l-  J')'>-fflLPHA*AR^fllN('T«J'i-flJNf 
r  TA~vv^™vxj\fy*wr^vxjw»xyy^V"'"vxjwy"  »r^vx'Xywyy'HwriiM~MLio\iyi^vxjwy*ij'*v 

1  "^A 
X  OO  m 

xx5*#*xyy»VwnvXfwy    ▼nvx*Xf*#yy    vw^vx    Xf»#y'vinivi    x#M*x/y*VTiivx  Xfwy'wfiv 

1  ^7 

cxyuy/     *^L.rnH*Moovw'^vx    xf^/wcrivx     xf*i«xyy^vwnivx    xyoy  wnivxfwyyy/v^.^iyc 

u  o  • 

J  I.  A  y 

X  5  ^  • 

FVYsTrVW/"?-  n4-VNfT-J-H'>l*fVNfT.  O-l-VfJfT  -  l-fll^*  ALPHA*  ftR^fVNfl-  I'l  +  VN 
f»i*vv"'»v  ifwy'Ti^v  x^w»x>y*x"i^vxfwy~»i^vxfv'xy/«Hu.rnH*MDovwrivxfwy' 

1  40 
X  "f  u  • 

lfT«J+lJ)*CVNCT.J)-VNCT«J+lJi— tVNtI.J-l)  +  VN(I.J)J*CVN(T.  J— 1  ^ -fYNC  T.J 

xvx9^     ^'^^^V*<^Kxy^y     *i~vx9w*xy.y    V"»'vxyw    xy'T'~vxywyy^v»i^vxjw  xy"w™vxpw 

141- 

X  ~  X  • 

2)l-ALPHA*AB<iCVNC  I.  J-1  )  +  VNCI«J)>*CVNCI.J-l^-VNCI«J))^/(4-»DELY") 

14  7- 

X  T  t  • 

VISX=  NU*(CUN(I+l«J)-2-»UN(I.J)+lINCI-l.J>)/nELX**2* 

wx**^       j^w^vvw^v^-'X^wy     *.»^wi^vxj*/y»i^i»vx  xjwyyrUL_4.rt^^4L.' 

143. 

1                      fUNCI*J+l>-2.*UNCI«J')  +  UNfI.J-l>'>/0ELY*«2') 

144- 

VISY=  NU*CCVNCI-i-i.J^-2-*VNtI«J>+VNCI-l.JJ3/0ELX**2-»' 

1  4<). 

1  CVNCr.J+l)-2-*VNCT.J)*VNCl.J-l))/DFLY**2) 

X                                 xv*^xX9w*xy     ^•^Tiivxywy'wi^vxyw  xyyriyi_wi^^*-y 

146  - 

U(I.J)=  UNCI.      ♦DELT*(  C  P(T  .  J  >- P(  I  ♦!  .  J">  ■)*  RDX   ♦   GXCI^-FUX-FUY*V  IS  X) 

147- 

1100 

VtI.JJ=  VNC I« J)+DELT*CCPCI« J)-PCI. J+1))*RDY  ♦   GYCI) -FVX-F VY*VIS Y) 

*xxy^^         wt^xxyv^'wiL.  w  i^xX'   XX  y  w  ^      I  xxyw'xyy^'^  w  ■             v'vxy      ■              •  v<*vx«jiy 

148. 

X  ~  7  • 

c  ♦  ♦ 

~  ~ 

SET   BOUNDARY  rONDlTTONS 

150. 

1  SI 

X  ^  X  • 

2000 

CDNTI NUF 

152- 

HNC1)=  HN(2) 

X  ^  J  • 

HNCINAX)=HNCIM1) 

X  3  ^  • 

JTC1)=JTC2) 

X  ^  7  • 

JT(I« AX)=JT(IM1) 

X  7  O  • 

C       LEFT  WALL  RIGID  AND  SLIP  FREE. 

1  57 

X  3  f  • 

C       RIGHT  WALL  CONTINUOUS  OUTFLOW. 

1  5  R  . 

X  3  O  • 

DO  2200  J^l.JMAX 

159- 

X  V  ^  • 

UC1»J)=0«0 

1  AO  . 

X  O  V  • 

VClf J>=VC2t J) 

X  ox  • 

IFCITER.GT.0>  GO  TO  2200 

1 

U<I«1,J)«U(I«2» J) 

1  A  ^- 

X  O  J  • 

VCIMAXf J)=VCIMlt J) 

1  A4 

2200 

CONTINUE 

1  A5  - 

C       TOP  MALL  CONTINUOUS  OUTFLOW. 

1  A  A  - 

X  o  o  • 

C       BOTTOM  WALL  RIGID  -  WITH  FRICTION 

1  A7  - 

X  O  f  • 

DO  2500  I=1,IMAX 

1  A  R  - 

X  w  o  • 

IFCITER.GT.O)  GO  TO  2400 

1  A9- 

X  O  7  • 

Vd.JMl  )=V(I,  JM2) 

1  70- 

X  1  V  • 

U(I,JMAX)=U(I,JN1) 

171- 
X  r  X  • 

2400 

VCI,1)=0.0 

17?- 
X  f  c  • 

2500 

U(I,l )=UClf2)*C1.0-2.0*FRCCI)) 

171- 

X   f  J9 

* 

174- 

X  f  ~  • 

C  *  ♦ 

FREE    SURFACE   BOUNDARY  CONDITIONS 

175 

X  1  J  • 

* 

1  76  - 

X  1  V  • 

DO  265J  I  =  2,I«^1 

1  77 
X  I  f  • 

JT1=JTC I) 

1  7  R 

X  f  o  • 

IFCJTCI  +  D.LT.JTd))  U(  I,  JT1)  =  UCI  .  JTl-1) 

179- 

X  1  7  • 

VCIf  JT1)=  VCI,JT1-1)-DELY*R0X*CUCI»JT1  )-UCI-l,  JTD) 

1  fiO 

X  o  u  • 

2650 

UCItJTl+l)=  UCIfJTli 

1  R  1 

X  O  X  • 

GO  TO  KRETt(3000t4100»4280) 

X  O  ^  • 

3000 

CONTINUE 

183. 

♦ 

18^. 

C  *  » 

HAS  CONVERGENCE  BEEN  REACHED 

185, 

♦ 

186. 

IFCFLG. EO.O.)  GO  TO  4000 

47 


187.  ITER=ITER+1 

188.  IFCITER.lt. 500)   GO  TO  3050 

189.  IFCCYCLE.LT.10>  GO  TO  4000 

190.  T=1.E+10 

191.  GO  TO  4000 

192.  3050  FLG=0.0 

193.  ♦ 

194.  C  *  ♦   COMPUTE  UPDATED  CELL  PRESSURE   AND  VELOCITIES 

195.  * 

196.  JBl=2 

197.  DO  3500   1  =  2  , IWl 

198.  JT1=JTCT) 

199.  00   3500  J=2fJTl 

200.  IFC JT l.EQ.JBl)   GO  TO  3060 

201.  IFCJ.NE.JBl    .AND.   J.NF.JTl)  GO   TO  3200 

202.  IFCJ.EQ.JTI)  GO  TO  3100 
^03.  GO  TO  3200 

204.  3060  CONTINUE 

20  5.  F=VCI  »J)*-DELT*RDX«CUCI,  J)-UCI-1,  J)) 

206.  DFDP=OE  LT*RDy*C 1.0  +  2.0 ♦DEL Y**2*«DX**2) 

207.  DELP1=-F/DFDP 

208.  3100  ETA=0ELY/CHNCI)-C'=L0ATCJT1)-2.5)*DELY) 

209.  DELP=C1.0-ETA)*PCIf JTi-l)-PCIf JT1> 

210.  IFC J8 l.EQ. JTl)  OELP=0.5*CDELP+OELP1) 

211.  GO  TO  3300 

212.  3200  0=RDX*(UC  I.  J)-UCI-1,  J))+R0Y*CVCI,  J)-VCI,  J-D) 

213.  IFCA3 SCD/DZRO).&E.EPSI>  FLG=1.0 

214.  OELP=  -BETA*0 

215.  3300  PCI tJ  )=PCIt J)+OELP 

216.  UCItJ)=UCIt J)*DELT*RDX*OELP 

217.  UCI-1 ,J >=UCI-lf J)-DELT*RDX*OELP 

218.  VCI, J>=VClf J>*OELT«ROY*DELP 

219.  3500  VCI tJ-l)=VC ItJ-l)-DELT»ROY*DELP 

220.  GO  TO  2000 

221.  4000  CONTINUE 

222.  » 

223.  C  »  »  COMPUTE  NEW  POSITION  FOR  TOP  SURFACE 

224.  ♦ 

225.  DO  4100  1=2, IMl 

226.  JT1=JTCI) 

227.  HV=  RDY*CHNCI)-FL0ATCJT1-2)«DELY) 

228.  UAV=   0.  5*CU(I-1,  JTl)  *  UCIfJTD) 

229.  HCI)=  HN(I)*FVOLl/FVOL+OELT«CHV*VCIf JTl )+C1.0-HV)*VCT,JTl 

230.  l-l)-0.5*ROX*CUAV*HNCI*l)*GA««A*ABS(UAV)*CHNCI)-HNtI*l)> 

231.  2  -UAV«HNCI-1)-GA««A*ABSCUAV)»CHNCI-1>-HNCI)))) 

232.  4100  CONTINUE 

233.  * 

234.  C  ♦  ♦  FIND   LEADING  AND  TRAILING  EDGES  OF   AVALANCHE  AND  LD .  EDGE  VELOCITY 

235.  ♦ 

236.  LDEG1=L0EG 

237.  I=I«2 

4400  IFCHCD.GT.DELW)  GO  TO  4500 

239.  1=1-1 

240.  GO  TO  4400 

241.  4300  LOEG'I 

242.  I«2 

243.  4600    IFCHC I) .GT.OELM)  60  TO  4700 

244.  1=1*1 

245.  GO  TO  4600 

246.  4700  KTE&=I 

247.  IFCLDEG.EO.LDEGl)   GO  TO  4800 

248.  IFCCYCLE.GT.O)  UEOG=DELX/TC 

249.  IFCCYCLE.EQ.O)  UEDG=5.0 

250.  TC«DELT 
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251.  INFL0=1 

252.  IFCUEDG.GT.UEOGl)  UEDG1=UEDG 

253.  00  ^750  I=KTEG,LDEG 

254.  4750  H(I)=H( I)*WCI-1)/WCI> 

255.  GO  TO  4910 

256.  4800  TC=TC+DELT 

257.  INFL0  =  IKFLO4-l 

258.  4910  CONTINUE 

259.  ♦ 

260.  C   ♦  ♦   CALCULfME   CELL   IN  WHICH  SURFACE   IS  LOCATEn   AND  UPDATE  ARRAY 

261.  * 

262.  DO   4250   1=2, IWi 

263.  IFCHt I) .LT.DELM)  HCI)=0.0 

264.  JT(  I)  =  INT(H(n*RDY+1.0E-6)*2 

265.  IF( JT(I).GT.J«1)  JT(I)=JH1 

266.  4250  CONTINUE 

267.  ASSIGN  4280   TO  KRET 

268.  GO  TO  200G 

269.  4280  CONTINUE 

270.  * 

271.  C  ♦   ♦  CALCULATE   TOTAL   ^LUID  VOLUME 

272.  ♦ 

273.  FVOL=0.0 

274.  00  4300  I=KTEG,LDEG 

275.  4300  FVOL=FV0L+H(I )»DELX*W(I) 

276.  IF( CYCLE. EQ.O)  FV0L1=FV0L 

277.  ♦ 

278.  C  *  ♦  ADVANCE   U,V,H  ARRAYS. 

279.  * 

280.  UPAX=VWAX=0.0 

281.  DO  4900  I=1,I«AX 

282.  DO  4900  J=1,JMAX 

283.  IFCABSCUd,  J)).GT.1.0E4^04)  U(I,J)  =  0.0 

284.  UNCI, J)=U(I  ,J> 

285.  TF( A3S(V(I, J)).&T.1.0E*04)  V(I,J)=0.0 

286.  VN( I, J)=VCI tJ) 

287.  IFCABSCPd,  J)).LT.i.QE-163  P(I,J)  =  0.0 

288.  490C  HNCI)=HCi; 

289.  DO  4950  I=KTEG,LDEG 

290.  DO  4950  J=2,JM1 

291.  UT=ABSCUN(I,J)) 

292.  VT=ABSt VN(I,J)) 

293.  IFC UT  .GT.UMAX)   UHAX  =  UT 

294.  4950   IF( VT.GT.VMAX)  VMAX=VT 

295.  * 

296.  C  ♦  »   LIST  VELOCITY,   PRESSURE,   AND  SURFACE  POSITION 
2  97.  * 

298.  5000       WRITE(108,49)   C YCLE , ITER  ,0 EL T , T , F VOL , U W AX , U£ OG ,LDEG 

299.  IFCCYCLE.NE .ICPRT)  GO  TO  6000 

300.  IF(CYCLE.EQ.ICPRT)  GO  TO  5030 

301.  5030   ICPRT=ICPRT  ^  INT(CWPRT) 

302.  5060  CONTINUE 

303.  WRITE(108,47) 

304.  DO  5250    1  =  1  ,IMAX 

305.  JT1=  JTCI) 

306.  JT2=JT1*1 

307.  00  5250  J=1,JT2 

308.  WRITE C108, 48)  I , J, UCI t J> t VC I i J) ( I t J) f HCI > • JTl 

309.  5250  CONTINUE 

310.  IFCIND.EQ.2)  GO  TO  6520 

311.  IFCIND.EQ.3>  GO   TO  6530 

312.  IFCIND.EQ.4)  GO  TO  6540 

313.  ♦ 

314.  C  ♦  ♦  RECOMPUTE  CONTROL  PARAMETERS. 
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315. 
316« 
317. 
318. 
319. 
320. 
321. 
322. 
323. 
324. 
325. 
326. 
327. 
328. 
329. 
330. 
331. 
332. 
333. 
334. 
335. 
336. 
337. 
338. 
339. 
340. 
341. 
342. 
343. 
344. 
345. 
346. 
347. 
348. 
349. 
350. 
351. 
352. 
353. 
354. 
355. 


oiFO: 
oiFo: 
OiFO: 
OlFO : 
OIFO : 
OiFo: 
OlFo: 
oiFO: 
OiFO: 
oiFo: 
OlFO: 
oiFo: 
oiFo: 
OlFO: 
OiFo: 
OlFO: 
oiFO: 
oiFO: 
OlFO: 
OiFO: 
oiFo: 


6000   IFCCYCLE.EQ.O)  GO  TO  6300 

DTX=DELX/U«AX 

DTY=OELY/VWAX 

DELT=A«INlC0TX,DTY)/3. 0 

IFCITER.lt. 10)  0ELT=1. 5*0ELT 

IFCNU-1 .E-6.lt. 0.0)  GO  TO  6300 

DET=C  DELX*DELY)»*2/C2.*NU»CDELX**2-«-0ELY*«2)) 

IF(DELT.LT.DET)  go  TO  6300 

DELT=0.9*DET 
6300  T=T*DELT 

IFCCYCLE.EQ.O)   GO  TO  6400 

DAX=U«AX*DELT/DELX 

DAY=VMAX*DELT/DELY 

ALPHA=1 .35* AWAXlCDAXtDAY) 

ifc alpha. gt. 1.0)  alpha=0.95 
Gamma^alpha 

8eta=0«g/c2.*delt*crdx*«24-rdy«*2)) 
c  *  ♦  test  for  program  termination. 

♦ 

6400  IFCT. GT.TWFIN)  IN0=2 

IFCHC I«2).GT.0ELM)  IN0=3 

IFCUEDG.lt. 0.05*UEDG1)  IN0=4 

IFCINFLO.EQ.50)  IND=4 

IFCIND.GT.l)  GO  TO  650U 

IF(CYCLE.LT.3)  GO   TO  6440 

AA=1.0+20.0*EXPC-1.25*UEOG) 

DO  64  30  I  =  2tim 
6430  FRCCI )=AA*FRC I) 
6440  CYCLE=CYCLE*1 

GO  TO  lOCO 
6500  T=r-DELT 

GO  TO   5  060 
6520  WRITECl08t82) 

GO  TO  6600 
6530  WRITEC108,83) 

GO  TO  6600 
6540  WRITEC108,84) 
6600  STOP 

END 
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Agriculture  —  CSU,  Fort  Collins 
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